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ABSTRACT 


The  interaction  of  the  northward  component  of  planetary  rotation  and  the  east- 
west  Reynolds  stress  affects  the  isotropy  of  the  integral  scale  turbulence  in  the  upper 
ocean  by  redistributing  turbulent  kinetic  energy  (TKE)  among  the  components.  This 
"rotation  stress"  mechanism  is  incorporated  into  a  vertically  integrated  model  of  the 
ocean  mixed  layer.  Simulations  of  Ocean  Weather  Stations  P  (50°N,  145°W)  and  N 
(30#N,  140° W)  are  used  to  compare  this  model  with  the  Garwood  (1977)  model  and 
with  observations.  The  significant  effect  is  the  augmentation  (for  easterly  winds)  or 
reduction  (for  westerly  winds)  of  the  ratio  of  vertical  to  horizontal  TKE.  The  rate  of 
entrainment  is  affected  by  the  change  in  the  vertical  convergence  of  TKE  at  the 
interface  between  the  mixed  layer  and  the  pycnocline.  Rotation  stress  significantly 
alters  the  mixing  on  diurnal  to  synoptic  time  scales  during  late  winter  and  early  spring. 
With  rotation  stress,  retreat  events  occur  more  frequently,  and  the  mixed  layer  depth 
change  during  retreat  is  10-30%  greater  than  without  rotation  stress.  Typically,  the 
ratio  of  vertical  to  total  TKE  is  three  times  larger  when  rotation  stress  is  included  and 
the  dissipation  enhancement  of  Garwood  (1977)  is  neglected.  The  resulting  TKE 
distribution  is  more  isotropic  and  in  better  agreement  with  laboratory  results  for 
neutrally  stratified  shear  flow's.  This  study  demonstrates  the  need  for  measurements  of 
the  TKE  budget  in  the  upper  ocean  to  confirm  these  findings  and  to  further  test  the 
hypotheses  of  TKE  models  in  oceanic  applications. 
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I.  INTRODUCTION 


This  dissertation  is  a  study  of  one  process  by  which  planetary  rotation  affects  the 
turbulent  kinetic  energy  (TKE)  budget  in  the  well-mixed  layer  of  the  upper  ocean.  The 
process  is  the  coupling  of  the  east-west  and  vertical  TKE  components  by  the 
interaction  of  the  northward  component  of  planetary  rotation  and  the  east-west 
Reynolds  stress.  This  results  in  a  redistribution  of  TKE  between  the  horizontal  and  the 
vertical  TKE  components. 

Turbulent  kinetic  energy  is  primarily  generated  near  the  sea  surface  by  the 
turbulent  fluxes  of  momentum  and  buoyancy.  The  entrainment  of  nonturbulent  fluid 
into  the  mixed  layer  occurs  at  the  base  of  the  mixed  layer.  Thus,  the  rate  of  mixed 
layer  deepening  depends  on  the  vertical  convergence  of  TKE  at  the  interface  between 
the  fully  turbulent  ocean  planetary  boundary  layer  (OPBL)  and  the  underlying 
pycnocline  (Garwood,  1977).  The  average  rate  of  vertical  turbulent  transport  through 
the  OPBL  to  the  interface  (entrainment  zone)  is  proportional  to  the  root  mean  square 
vertical  turbulent  velocity  averaged  over  the  mixed  layer.  Both  the  magnitude  of  the 
vertical  TKE  averaged  over  the  mixed  layer  and  the  magnitude  of  the  total  TKE  are 
important  to  the  dynamics  of  the  upper  ocean  mixed  layer. 

Through  the  rotation  stress,  that  is  the  interaction  of  the  northward  component 
of  planetary  rotation  and  the  east-west  Reynolds  stress,  the  mixed  layer  dynamics  is 
predicted  to  depend  on  latitude  and  wind  direction.  The  transfer  of  energy  increases 
the  east-west  component  of  TKE  and  decreases  the  vertical  TKE  for  westerly  winds 
(eastward  surface  Reynolds  stress).  The  decrease  in  the  vertical  TKE  reduces  the 
convergence  of  TKE  to  the  entrainment  zone  and  reduces  the  rate  of  mixed  layer 
deepening.  Conversely,  the  rate  of  deepening  is  augmented  for  easterly  winds 
(westward  surface  Reynolds  stress). 

A.  STATEMENT  OF  THE  PROBLEM  AND  PLAN  OF  ATTACK 

This  study  will  explore  the  diurnal,  synoptic  and  annual  variations  in  mixed  layer 
dynamics  which  result  from  rotation  stress  for  mid  latitude  mixed  layers  in  both 
easterly  and  westerly  wind  regimes.  A  vertically  integrated,  second-order  closure  model 
of  the  OPBL  will  be  used  to  calculate  the  time-dependent  mixed  layer  depth  and 
temperature  and  the  distribution  of  energy  between  the  horizontal  and  vertical  TKE 
components. 


The  time  scales  for  rotation  stress  are  the  earth's  rotation  time  scale  and  the  time 
scale  for  variations  in  the  wind  direction.  The  diurnal  to  synoptic  time  frame 
encompasses  these  time  scales.  On  diurnal  time  scales,  the  mixed  layer  frequently 
retreats  in  response  to  the  solar  insolation  then  deepens  at  night.  On  synoptic  time 
scales,  increases  or  decreases  in  surface  wind  stress  and  solar  insolation  caused  by  the 
passage  of  atmospheric  low  pressure  systems  can  cause  the  mixed  layer  to  deepen  or 
retreat.  Following  the  passage  of  the  disturbance,  the  mixed  layer  will  readjust  to  the 
local  atmospheric  conditions  and  may  return  to  approximately  the  same  depth  which 
existed  prior  to  the  event.  For  these  time  scales,  the  rotation  stress  may  significantly 
alter  the  frequency  and  vertical  scale  of  these  shallowing  and  deepening  events 
depending  on  the  variations  in  wind  direction  throughout  the  event.  This  is  particularly 
likely  for  cases  when  the  wind  shear  production  and  the  surface  buoyant  damping  of 
TKE  are  nearly  equal  and  the  mixed  layer  is  deep  0<100m). 

Previous  studies  of  the  mid  latitude  OPBL  have  demonstrated  that  the  annual 
cycle  of  mixed  layer  shallowing  and  deepening  is  poorly  simulated  unless  the 
dissipation  is  increased  during  times  of  increased  buoyant  production  of  TKE  (see 
Zilitinkevich,  1979  and  Garwood,  1979  for  reviews).  Consequently,  the  dissipation  in 
mixed  layer  models  has  been  augmented  in  various  ways  to  prevent  excessive 
wintertime  mixed  layer  deepening.  The  need  for  this  enhancement  of  the  dissipation 
indicates  that  physical  processes  that  are  important  for  mixed  layer  dynamics  are  not 
being  included  in  the  models  or  that  the  processes  are  being  modeled  incompletely. 
Although  rotation  stress  has  been  shown  to  be  important  for  the  equilibrium  mixed 
layer  depth  in  the  tropics  (Garwood  el  al.,  1985b),  it  has  not  previously  been  included 
in  a  model  of  the  time-dependent  mixed  layer.  For  westerly  wind  regimes,  it  is  possible 
that  the  reduction  in  the  entrainment  rate  due  to  rotation  stress  may  reduce  or 
eliminate  the  need  for  enhanced  dissipation.  Conversely,  it  may  not  be  possible  to 
achieve  a  neutral  equilibrium  depth  for  the  OPBL  based  on  local  surface  forcing  for 
regions  of  predominately  easterly  winds.  In  the  latter  case,  a  net  downward  heat  flux 
or  a  three-dimensional  circulation  may  be  required  to  achieve  a  closed  annual  c\  cle  for 
the  mixed  layer  depth  in  easterly  wind  regimes. 

The  next  section  of  this  dissertation  gives  the  background  necessary  for  this 
study.  The  governing  equations  are  given,  and  the  terms  are  defined  and  discussed. 
Then  a  physical  mechanism  is  proposed  to  explain  the  rotation  stress.  The  turbulence 
closure  problem  and  the  resulting  need  for  parameterizations  arc  addressed.  The 


parameterizations  for  rotation  stress  and  dissipation  are  discussed  briefly.  The 
parameterization  of  the  pressure-rate  of  strain  term  is  discussed  in  some  detail  because 
it  contains  terms  which,  in  theory,  could  counter  the  effects  of  the  rotation  stress 
terms.  In  the  final  section,  the  hypotheses  and  tests  of  these  hypotheses  are  put 
forward. 

Chapter  II  contains  the  literature  review  in  three  sections.  The  first  section  deals 
with  the  development  of  vertically  integrated  (bulk)  mixed  layer  models.  The  second 
section  reviews  previous  work  concerning  the  effects  of  rotation  on  turbulent  boundary 
layers.  The  effects  of  streamline  curvature  are  also  reviewed  since  they  are 
mathematically  identical  to  those  for  rotation.  Finally,  the  problems  involved  in 
parameterizing  the  pressure-rate  of  strain  terms  are  reviewed  since  this 
parameterization  affects  the  results  of  this  work. 

In  Chapter  III,  the  vertically-integrated  equations  for  the  TKE  components  are 
developed  and  discussed.  The  nondimensional  equations  are  derived  for  the  models 
that  will  be  used  in  this  study.  The  TKE  budget  is  examined  to  highlight  the 
differences  between  the  processes  that  are  explored  in  this  dissertation.  The  models  are 
compared  using  the  nondimensional  framework  introduced  by  Garwood  (1977)  and 
used  for  model  intercomparison  in  Garwood  (1979). 

Chapter  IV  contains  a  description  of  the  data  sets  extracted  from  Ocean  Weather 
Station  (OWS)  data  archives,  which  will  be  used  for  the  simulations  in  this  study.  The 
results  of  the  annual  simulations  are  discussed  for  synoptic,  seasonal  and  annual  time 
scales.  The  effects  of  adjusting  the  tuning  constants  are  also  explored.  Chapter  V 
contains  a  summary  of  the  conclusions  and  recommendations  resulting  from  this  study. 

B.  BACKGROUND 

1.  Governing  equations  and  definitions  of  terms 

The  upper  ocean  mixed  layer  is  approximately  horizontally  homogeneous  if 
the  smallest  term  involving  vertical  derivatives  is  sufficiently  larger  than  the  largest 
term  involving  horizontal  derivatives.  To  formally  determine  if  this  criterion  is  met,  all 
the  terms  in  the  TKE  budget  must  be  approximated  systematically  and  the  order  of 
approximation  at  which  different  terms  are  eliminated  must  be  analyzed.  Pedloskv 
(1979)  applies  this  technique  in  detail  to  derive  the  shallow  water  equations.  However, 
Businger  (1982)  points  out  that  the  assumption  of  horizontal  homogeneity  is  usually 
valid  if  h'L  <  <  1 ,  where  h  is  a  length  scale  for  the  vertical  and  L  is  a  horizontal  length 
scale.  Since  the  TKE  in  the  upper  ocean  responds  to  the  local  atmospheric  forcing, 


one  relevant  horizontal  length  scale  in  the  mid  latitudes  would  be  that  of  extra  tropical 
cyclones.  This  length  scale  is  Vg’ii/i  —  1000  km,  where  g'  =  g(po-p)  po  is  the 
reduced  gravity.  This  length  scale  is  the  baroclinic  Rossby  radius  for  the  atmospheric 
synoptic  scale  with  a  vertical  scale  height  Z~  10km  (Pcdlosky,  1979).  The  mean 
dynamical  variables  in  the  ocean  have  a  synoptic  spatial  scale  of  Lf  =  Vg  Z0”  f~  30km 
where  ZQ  is  the  depth  scale  of  the  permanent  pycnocline  0(500m).  For  L  equal  to  Lf 
or  Ls,  h'L<  <1  and  horizontal  homogeneity  is  a  valid  assumption.  Since  L$>Lr, 
horizontal  homogeneity  is  probably  a  better  assumption  for  the  TKE  budget  than  for 
the  mean  flow. 

A  set  of  three  equations  describes  the  TKE  budget  in  the  OPBL.  The 
equations  are  derived  using  Reynolds  decomposition  (Tennekes  and  Lumley,  1972). 
The  origin  of  the  coordinate  system  is  at  the  sea  surface,  x  is  eastward,  y  is  northward 
and  z  is  upward.  For  the  components  of  TKE,  which  are  the  normal  Reynolds 
stresses,  the  equations  are: 
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Here,  the  lower  case  letters  represent  fluctuations  of  the  variables  about  the  means, 
which  are  shown  as  upper  case  letters.  The  overbars  represent  ensemble  averages. 

Turbulence  is  isotropic  if  the  statistically  averaged  properties  of  the  flow  are 
invariant  under  arbitrary  rotation  and  reflection  of  the  coordinate  system  (Hinze, 
1975).  A  consequence  of  isotropy  is  that  the  TKE  would  be  distributed  equally  among 
the  components  (u2,  v2  and  w2)  in  (1.1),  (1.2)  and  (1.3).  The  OPBL  is  inherently 
anisotropic,  as  are  all  turbulent  shear  flows,  since  the  mean  shear  imposes  a  preferred 
direction.  The  generation,  transport  and  redistribution  of  turbulence  in  such  mixing 
layers  is  a  direct  consequence  of  this  anisotropy.  The  TKE  distribution  among  u2,  v^ 
and  w'  is  anisotropic  since  the  the  shear  production  (G.l  and  G.2)  generates  TKE  only 


in  the  u2  and  v-  components  (1.1  and  1.2),  and  the  buoyancy  production  and  damping 
(H)  generates  and  damps  TKE  only  in  the  vertical  w2  component  (1.3).  At  the  scales 
of  the  production  terms,  the  TKE  components  associated  with  the  shear  production 
can  be  a  factor  of  two  larger  than  the  other  components  (Townsend,  1976).  The 
processes  represented  by  the  terms  in  (1.1  -  1.3)  can  be  categorized  with  respect  to  the 
anisotropy  of  the  TKE  distribution.  A  source  of  TKE  for  the  OPBL  is  the  shear 
production  due  to  the  shear  of  the  mean  current  and  the  surface  wind  stress,  terms  G.  1 
-  G.2.  Also  included  in  these  terms  is  a  sink  of  TKE,  the  entrainment  stress  which  is 
due  to  the  entrainment  of  nontur’oulent  fluid  at  the  base  of  the  mixed  layer.  The 
damping  due  to  the  entrainment  shear  reduces  only  the  horizontal  TKE.  The 
buoyancy  flux,  H.  represents  both  the  surface  buoyancy  flux,  which  can  be  a  source  or 
sink  of  vertical  TKE  and  the  entrainment  buoyancy  flux,  which  is  a  sink  of  vertical 
TKE.  This  production  and  damping  by  buoyancy  flux  increases  the  anisotropy 
because  it  alters  only  the  vertical  component  of  TKE  (1.3). 

The  convergence  of  TKE  (P.l  -  P.3)  at  the  interface  between  the  turbulent 
OPBL  and  the  dynamically  stable  pycnocline  is  a  major  determinant  of  the  entrainment 
rate  (Garwood,  1977).  The  entrainment  rate  is  the  rate  at  which  nonturbulent  water  is 
entrained  into  the  OPBL,  and  hence  the  rate  at  which  the  mixed  layer  deepens.  The 
convergence  of  TKE  depends  on  the  vertical  turbulent  velocity,  which  scales  with  the 
square  root  of  the  vertical  TKE.  Thus,  the  anisotropy  of  the  TKE  distribution  and  the 
magnitude  of  the  TKE  are  both  important  to  mixed  layer  dynamics.  Processes  that 
affect  the  isotropy  of  the  TKE  distribution  by  transferring  TKE  between  the  horizontal 
and  vertical  components  can  significantly  impact  mixed  layer  dynamics  even  though 
they  do  not  directly  alter  the  total  TKE.  The  terms  that  represent  such  processes  are 
n.l  -  11.3,  and  R.l  and  R.2. 

The  interaction  of  planetary  rotation  and  the  Reynolds  stresses  is  represented 
by  (R.l  -  R.4).  These  rotation  stresses  depend  on  latitude  and  wind  direction.  The 
northward  (lly)  and  vertical  (ft 2)  components  of  planetary  rotation  are: 

£ly  =  Qcos6, 

£  =  2f2sin6, 

where  0  is  the  latitude.  The  surface  Reynolds  stresses  (uw(0)  and  vw(0))  are  the 
vertical  flux  of  momentum  due  to  the  wind  stress, 

■f  =  -  p(uwftO)  +  vw(  0)  fcy)  . 


(1.4) 


The  rotation  stress  is  the  work  done  on  the  turbulent  flow  by  the  planetary  rotation. 
For  westerly  surface  winds  (eastward  or  positive  tx).  the  rotation  stress  transfers  TKE 
from  the  vertical  to  the  horizontal  components.  This  reduces  the  vertical  TKE  and 
reduces  the  vertical  transport  of  TKE  which  decreases  the  rate  at  which  the  mixed  layer 
deepens.  If  the  primary  turbulence  generation  mechanism  is  wind  shear  production 
(forced  convection),  the  anisotropy  of  the  turbulence  would  be  increased,  since  both 
the  primary  production  (G.l  and  G.2)  and  the  rotation  stress  (R.l  and  R.2)  would  be 
increasing  the  energy  in  the  horizontal  component.  For  forced  convection  regimes  with 
easterly  winds,  the  entrainment  rate  is  enhanced  and  the  turbulence  becomes  more 
isotropic  through  the  rotation  stress  interaction. 

The  pressure-rate  of  strain  (II.  I  -  II. 3)  is  the  cross-correlation  of  the  dynamic 
pressure  fluctuations  and  the  turbulent  velocities.  The  dynamic  pressure  fluctuations 
result  from  local  and  remote  variations  in  the  mean  or  turbulent  velocity  field.  Any 
convergence  or  divergence  of  velocity  can  set  up  an  adverse  gradient  of  fluctuating 
pressure.  The  pressure-rate  of  strain  is  a  reaction  of  the  turbulent  flow’  to  variations  in 
the  turbulence  and  to  variations  in  the  mean  flow.  The  pressure-rate  of  strain  due  to 
the  turbulence-turbulence  interaction  is  comprised  of  the  turbulent  velocity  gradient 
and  the  dynamic  pressure  fluctuations  that  result  from  variations  in  the  turbulent 
velocity.  The  primary  effect  of  this  pressure-rate  of  strain  interaction  is  to  reduce  the 
anisotropy  of  the  turbulence  by  causing  a  transfer  of  TKE  from  the  component  in 
which  most  of  the  production  is  taking  place  to  the  other  TKE  components 
(Townsend.  1976).  For  forced  convection  regimes  in  w’hich  shear  production  (G.l  and 
G.2)  is  the  primary  production  mechanism  for  the  OPBL,  this  implies  that  the  effect  of 
the  turbulence-turbulence  pressure-rate  of  strain  is  to  reduce  the  energy  in  the 
horizontal  components  and  increase  the  energy  in  the  vertical  component.  The  other 
part  of  the  pressure-rate  of  strain  is  the  turbulence-mean  flow  interaction,  which  is 
comprised  of  the  turbulent  velocity  gradient  and  the  dynamic  pressure  fluctuations  that 
result  from  variations  in  the  mean  velocity.  The  nature  and  effects  of  the  turbulence- 
mean  flow  part  of  the  pressure-rate  of  strain  are  dependent  on  the  nature  of  the 
particular  flow  regime.  However,  there  is  evidence  that  it  may  be  negligible  for 
geophysical  turbulent  boundary  layers  (see  Section  3). 

The  rotation  stress  (R.l  -  R.4)  and  the  pressure-rate  of  strain  (II.  1  -  11.3)  act 
to  redistribute  TKE  among  the  components  without  directly  influencing  the  total  TKE. 
This  can  be  seen  by  summing  (1.1),  (1.2)  and  (1.3)  to  obtain  an  equation  for  twice  the 
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total  TKE,  e  =>  u2  ^  v2  +  w2.  In  the  resulting  equation,  the  terms  R.  1  -  R.4  and  n.l  - 

11.3  sum  to  zero  so  the  total  TKE  is  not  affected  directly.  Because  many  ocean  mixed 

layer  models  use  only  the  total  TKE  equation  (see  Ziiitinkevich  ei  al.,  1979  and 

Garwood,  1979  for  reviews),  the  rotation  stress  and  pressure-rate  of  strain  terms  are 

not  included  in  those  models.  In  other  cases  (e.g.  Wyngaard,  19S1).  2ft  uw  is 

neglected  based  on  scaling  the  Reynolds  equations  rather  than  the  TKE  equations.  In 

the  Reynolds  equations,  the  equivalent  term  is  2ft  w.  This  scales  as 

_  fu*  h 

2CI  w  ~  —  , 

where  the  ratio  of  w  to  u*  is  derived  from  the  continuity  equation  and  u*  is  a 
characteristic  horizontal  turbulent  velocity  scale  based  on  the  surface  Reynolds  stress 
(1.4) 

u*2  *  t(0)/p~  0.01  m  s"2. 

With  the  boundary  layer  approximation,  h/Lf<  <1,  this  term  is  much  smaller  than  the 
other  Coriolis  term  fu,  since  w<  <u*  from  the  continuity  equation.  Therefore  it  is 
negligible  in  the  turbulent  momentum  equation.  However,  uw  and  vw  are  normally 
larger  than  uv  since  they  are  the  Reynolds  stresses  associated  with  the  surface  wind 
stress  (1.4).  Thus,  R.l  and  R.2  are  at  least  as  large  as  R.3  and  R.4  and  these  terms 
should  be  retained  during  the  Reynolds  averaging  process.  The  importance  of  the 
rotation  stress  terms  for  the  TKE  budget  can  be  determined  by  comparing  the  rotation 
stress  to  the  pressure-rate  of  strain  since  these  are  the  terms  that  causes  a  transfer 
energy  among  the  TKE  components  (1. 1-1.3).  The  ratio  of  the  pressure-rate  of  strain 
to  the  rotation  stress  is 
w(3p 

vdz  u*3/h _  ,,  _ 
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So  the  rotation  stress  scales  with  a  boundary  layer  Rossby  number  (Roh)  in  the  TKE 
equations.  For  u*“0.01  m  s’1,  which  corresponds  to  a  surface  stress  of  0.1  N  m'2  (or 
wind  speed  of  8  m  s'1  at  10  m),  and  h*  100  m.  Roh-  1.  Thus,  the  redistribution  of 
TKE  by  rotation  stress  can  not  be  neglected  based  on  scaling  arguments.  The  relative 
magnitude  of  the  rotation  stress  increases  with  decreasing  Rossby  number.  Therefore, 
rotation  stress  becomes  relatively  more  important  for  mixed  layer  dynamics  for  deep 
mixed  layers  and  low  wind  speeds. 

The  viscous  dissipation  (£*  D.l  +  D.2  +  D.3)  is 


where  i-  1,2,3  and  j»  1,2,3,  corresponding  to  the  x.  y,  z  axis  respectively.  The 
dissipation  is  assumed  to  be  isotropic,  so  that  energy  is  removed  equally  from  all  three 
TKE  components.  Approximately  80  to  90%  or  more  of  the  TKE  produced  by  G.l, 
G.2  and  H  is  dissipated  locally  by  D.l  -  D.3  (Tennekes  and  Lumley,  1972).  The 
unsteady  terms,  T.l  -  T.3,  may  augment  or  reduce  the  isotropy  depending  on  the  initial 
state  of  the  turbulence  and  on  the  signs  of  the  terms.  However,  these  terms  are  usually 
negligible  in  geophysical  planetary  boundary  layers  except  for  time  scales  shorter  than 
100  seconds  (deSzoeke  and  Rhines,  1976  and  Garwood  and  Yun,  1979  ). 

2.  A  Physical  Mechanism 

Turbulence  seems  to  exist  as  a  collection  of  interspersed  and  interlinked 
vortices  (eddies)  of  various  sizes  and  shapes.  Here,  the  term  eddy  means  a  spatially 
limited  distribution  of  vorticity  in  a  fairly  simple  pattern,  such  as  a  cylinder  or  roll  cell, 
or  a  sphere  (Townsend,  1976).  The  sizes  range  from  the  integral  scale,  which  is  the 
largest  scale  of  the  flow,  to  the  Kolmogorov  scale,  which  is  the  scale  of  the  dissipation 
of  TKE  by  molecular  viscosity.  The  rotation  stress  will  be  most  effective  at  the  largest 
scales  of  the  turbulent  flow  since  the  ratio  of  the  pressure-rate  of  strain  to  the  rotation 
stress  is  proportional  to  1  h  (1.5).  Thus  the  effect  of  rotation  stress  will  be  examined 
for  a  roll  vortex  with  the  rotation  axis  oriented  north-south  (along  the  y  axis  of  the 
reference  frame)  and  with  a  scale  normal  to  the  rotation  axis  that  is  comparable  to  the 
mixed  layer  depth,  which  is  the  integral  scale  for  the  OPBL  (Figure  1.1). 

Rotation  stress  is  the  w'ork  done  on  these  eddies  by  the  Coriolis  forces.  In  the 
Reynolds  stress  equations,  the  terms  which  represent  the  northward  component  of  the 
Coriolis  force  are: 
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dt 
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These  terms  represent  inertial  motion  in  the  x-z  plane  that  is  completely  analogous  to 
the  more  familiar  inertial  motion  in  the  x-y  plane.  For  westerly  winds  (du,dz>0),  this 
inertial  motion  induces  a  positive,  or  clockwise  (looking  in  the  positive  y  direction  or 
north),  circulation  or  vorticity.  Changes  in  vorticity  are  accompanied  by  changes  in 


size  of  the  vortex  (Tennekes  and  Lumley,  1972).  The  vertical  scale  of  the  roll  vortex  is 
reduced  for  westerly  winds  as  the  vorticity  is  increased.  This  reduces  the  vertical  scale 
of  mixing  (Figure  l.l).  For  easterly  winds  (<Ju/dz<0),  the  inertial  motion  induces  a 
negative  circulation  or  vorticity.  Thus,  the  vorticity  is  reduced  and  the  vertical  scale  of 
the  roll  vortex,  and  hence  the  mixed  layer  depth,  is  increased.  The  vorticity,  the 
distribution  of  TKE,  the  rate  of  mixed  layer  deepening  and  the  equilibrium  depth  of  the 
mixed  layer  may  all  be  altered  by  rotation  stress. 

Do  eddies,  roll  vortices  or  other  coherent  structures  actually  exist  in  the 
OPBL?  There  is  evidence  that  regular,  large-scale,  coherent  structures  exist  for  many 
turbulent  flows.  Intermittent  turbulent  activity  observed  in  the  atmospheric  boundary 
layer  is  suggestive  of  large-scale,  verticaliy-coherent  structures  (Shaw  and  Businger, 
1985).  Various  large-scale  structures  have  been  observed  in  laboratory  flows  (Cantwell, 
19S1).  These  include  the  Taylor-Goertler  vortex  which  is  a  well-documented  secondary 
circulation  that  is  found  in  rotating  laminar  and  turbulent  flows  (Johnston,  1972;  Hunt 
and  Joubevt,  1979).  Taylor-Goertler  vortices  are  cross-stream  roll  vortices  that  result 
from  the  destabilizing  effects  of  rotation  on  the  mean  flow.  Langmuir  circulations  have 
been  observed  in  the  upper  ocean  (see  review  by  Leibovich,  1983)  and  longitudinal  roll 
vortices  have  been  observed  in  the  atmosphere  (LeMone,  1973).  It  is  possible  that 
both  types  of  Langmuir  circulations  are  planetary  boundary  layer  scale  Taylor-Goertler 
vortices. 

Although  the  existence  of  these  large-scale,  evanescent  structures  is  becoming 
undeniable,  it  has  not  been  demonstrated  that  explicitly  simulating  these  structures 
improves  model  performance  (Saffman,  1980).  The  models  developed  and  used  in  this 
research  belong  to  a  class  known  as  bulk  or  vertically-integrated  models.  Thus,  the 
structure  of  the  turbulence  is  not  explicitly  accounted  for  in  these  models.  Models  of 
this  type  take  advantage  of  the  geometric  similarity  of  the  upper  ocean  structure  to 
simplify  the  mean  momentum  and  the  TKE  equations.  The  similarity  solution  used  in 
this  study  defines  the  upper  ocean  as  a  fully  turbulent,  and  therefore  well-mixed,  layer 
of  depth  h.  Thus  the  OPBL  is  homogeneous  in  mean  velocity  and  density.  This  implies 
that  the  Reynolds  stresses  are  constant  or  linear  in  z  (see  Chapter  III).  The  layer  is 
bounded  above  by  a  wind-wave  zone  of  depth  6W  <  <  h  and  at  the  bottom  by  a  very 
sharp  density  gradient  (the  entrainment  zone)  of  depth  8<  <h,  (Figure  1.2).  The 
similarity  structure  is  essentially  that  of  Niiler  (1975)  and  Garwood  (1977). 
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3.  The  Closure  Problem  and  Parameterizations 


The  equations,  (1.1)  -  (1.3),  for  the  TKE  components  in  the  OPBL  are  the 
equations  for  the  second-order  autocorrelations  of  the  turbulent  velocity  components. 
These  three  equations  contain  14  second-  and  third-order  correlations  that  are 
unknown.  These  include  the  shearing  Reynolds  stresses  (£Tuj),  the  pressure- velocity 
correlations  (puj),  the  vertical  transports  (wu;2  and  wp),  the  buoyancy  flux  (bw)  and 
the  viscous  dissipation  terms  (1.6).  The  inclusion  of  equations  for  these  unknown 
correlations  incorporates  more  unknown  correlations  of  equal  and  higher  order.  This 
problem  is  inherent  in  the  process  of  Reynolds  decomposition  and  constitutes  the 
turbulence  closure  problem  (Tennekes  and  Lumley,  1972).  The  solution  to  this 
problem  is  to  choose  a  particular  set  of  equations  and  "parameterize"  the  unknown 
correlations  in  terms  of  known  correlations  of  the  same  or  lower  order,  mean  variables 
and  boundary  conditions.  These  parameterizations  are  based  on  laboratory  and  field 
observations,  dimensional  analysis  and  physical  intuition.  The  order  of  the  correlations 
that  are  calculated  explicitly,  the  number  of  equations  that  are  solved  and  the 
parameterizations  of  the  terms  involving  unknown  correlations  varies  among  models  of 
planetary  boundary  layers  (Wyngaard,  1982). 

Some  of  the  models  for  the  OPBL  will  be  reviewed  in  Chapter  II.  The  models 
used  in  this  research  are  bulk  (vertically-integrated),  second-order  closure  models  as  in 
Garwood  (1977).  The  details  of  the  particular  closure  method  and  parameterization 
schemes  used  in  this  dissertation  will  be  discussed  in  Chapter  III.  However,  three 
processes  and  their  parameterizations  that  will  be  particularly  significant  for  this 
research  need  to  be  introduced  before  the  working  hypotheses  and  tests  of  the 
hypotheses  are  discussed.  These  are  the  dissipation,  pressure-rate  of  strain  and 
rotation  stress. 

This  research  is  concerned  with  the  rotation  stress  terms  (R.l  and  R.2)  that 
alter  the  ratio  of  vertical  to  horizontal  TKE  in  (1.1)  and  (1.3).  These  terms  affect  the 
distribution  of  TKE  between  the  east-west  and  vertical  TKE  components.  The  average 
vertical  turbulent  velocity  is  changed  which  affects  the  rate  of  mixed  layer  deepening  by 
altering  the  convergence  of  TKE  at  the  entrainment  zone.  The  parameterization  of 
these  terms  requires  no  assumptions  beyond  those  of  the  similarity  solution.  Garwood 
et  al.  (1985a)  derived  a  bulk  model  that  included  rotation  stress  for  a  non-entraining 
mixed  layer.  They  showed  that  the  equilibrium  depth  that  is  achieved  for  a  surface 
wind  stress,  T  ( 1 .4),  and  a  positive  downward  surface  buoyancy  flux  is 
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( i  +  c2a>) 

where  the  Obukhov  length  scale, 

r  V2 

L  =  - - -  ,  (1.8) 

0  'Ew(0) 

is  the  equilibrium  depth  for  a  balance  between  the  surface  stress  and  surface  buoyancy 
flux.  The  ratio  (<l>)  of  the  rotation  stress  to  the  surface  buoyancy  flux  is 

n  t 

a>  =  —u.  .  (i.9) 

pBw( 0) 

The  equilibrium  depth,  Lq,  reflects  the  vertical  elongation  or  contraction  of  mixed 
layer  vortices  in  the  x-z  plane  by  planetary  rotation  shown  in  Figure  1.1.  Rotation 
stress  effectively  modulates  the  damping  due  to  the  surface  buoyancy  flux  such  that  the 
equilibrium  depth,  Lq,  will  be  greater  for  easterly  winds  (westward  surface  stress)  or 
less  for  westerly  winds  (eastward  surface  stress)  than  Lo,  the  Obukhov  depth  (Figure 
1.3).  We  have  found  no  physical  solution  for  Lq  is  possible  for  <P<-l/2  because  total 
buoyant  damping  can  not  exceed  dissipation.  As  the  minimum  value  of  <J>  =  -l/2  is 
approached,  more  of  the  TKE  is  damped  by  surface  buoyancy  flux,  and  the  equilibrium 
depth  approaches  the  limit  Lq  =  7Lq.  For  positive  <I>,  westerly  winds,  Lq<  Lq  due  to 
the  increased  transfer  of  vertical  to  horizontal  TKE  by  rotation  stress.  The  increased 
horizontal  TKE  is  dissipated,  and  the  mixing  is  reduced.  Figure  1.3  shows  that 
Lq<O  =  0  .58)  is  half  the  value  of  Lq(<J>  =  0).  Thus,  for  westerly  winds  rotation  stress 
causes  an  enhanced  dissipation  of  TKE  and  results  in  a  reduced  equilibrium  depth. 

Planetary  rotation  has  previously  been  incorporated  into  a  bulk  model  of 
turbulence  in  the  OPBL  through  the  dissipation  time  scale  (Garwood,  1977).  If  the 
dissipation  rate  is  governed  by  the  rate  at  which  energy  is  supplied  by  the  large  scale 
turbulence,  then  dimensional  arguments  yield  e^e^'^h  (Tennekes  and  Lumley,  1972). 
This  assumes  that  the  time  scale  for  the  energy  transfer  is  the  "eddy  turnover  scale" 
(h.'Ve),  where  the  mixed  layer  depth  is  an  appropriate  length  scale  for  the  large-scale 
flow.  This  parameterization  has  been  shown  to  be  inadequate  for  slab  models  of  the 
OPBL,  since  it  allows  a  too  deep  winter  mixed  layer  in  annual  simulations  (see  reviews 
by  Garwood,  1979,  and  Zilitinkevich  et  al.,  1979).  Garwood  (1977)  argued  that 
planetary  rotation  provides  an  additional  time  scale  (I  f)  for  dissipation.  This  becomes 
important  for  deep  mixed  layers  for  which  I  fsh  v^  The  dissipation 
parameterization  becomes  £~e'(v/^"h+0-  This  enhancement  to  the  dissipation 
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parameterization  is  a  rather  ad  hoc  method  of  incorporating  planetary  rotation  into  a 
mixed  layer  model.  It  is  not  clear  what  physical  processes  are  responsible  for  the 
additional  dissipation.  For  westerly  winds,  rotation  stress  reduces  the  vertical 
convergence  of  TKE  into  the  entrainment  zone.  This  reduces  the  rate  of  mixed  layer 
deepening.  Thus,  rotation  stress  may  reduce  or  eliminate  the  need  for  the  planetary 
rotation  "enhancement"  to  the  dissipation  for  westerly  wind  regimes.  Furthermore, 
rotation  stress  represents  a  specific  physical  mechanism  rather  than  a  time  scale  based 
on  dimensional  arguments.  One  aspect  of  this  research  will  be  to  compare  the  effects 
of  rotation  stress  in  westerly  wind  regimes  to  the  rotational  dissipation  enhancement  of 
Garwood  (1977). 

The  pressure-rate  of  strain  (11.1  -  n.3)  is  the  only  process  other  than  rotation 
stress  that  redistributes  energy  among  the  TKE  components.  These  are  the  only  terms 
that  may  directly  counter  or  reduce  the  effects  of  rotation  stress.  The  expression  for 
the  fluctuating  pressure,  which  is  the  basis  for  parameterizations  of  the  pressure-rate  of 
strain,  contains  terms  that  are  similar  to  the  rotation  stress  terms  (R.l  -  R.4)  in  (1.1  - 
1.3),  but  with  opposite  sign.  These  terms  in  the  pressure-rate  of  strain 
parameterization  may  reduce  the  effects  of  the  explicit  rotation  stress  terms. 

The  parameterization  of  the  pressure-rate  of  strain  requires  an  understanding 
of  both  the  fluctuating  pressure  and  the  gradient  of  the  fluctuating  velocity.  As 
mentioned  above,  the  fluctuating  pressure  can  react  to  local  and  remote  variations  of 
the  turbulent  and  mean  flows.  However,  the  exact  effects  of  the  various  processes  that 
contribute  to  the  pressure  fluctuations  are  not  well  documented  by  experiment  nor  are 
they  well  understood  theoretically.  Thus,  parameterizing  the  pressure-rate  of  strain 
remains  an  active  research  area,  although  this  problem  has  been  apparent  for  more 
than  40  years.  It  is  beyond  the  scope  of  this  dissertation  to  address  the  problem  of 
parameterizing  the  pressure-rate  of  strain  in  general.  However,  a  realistic 
parameterization  of  the  term  is  needed  to  complete  the  model.  The  specific 
parameterization  of  the  pressure-rate  of  strain  that  is  used  for  this  research  will  be 
discussed  in  Chapter  III. 

The  possibility  that  the  pressure-rate  of  strain  may  tend  to  cancel  the  rotation 
stress  will  be  addressed  here.  An  expression  for  the  fluctuating  pressure  can  be 
obtained  by  taking  the  divergence  of  the  fluctuating  velocity  equations  and  making  use 
of  the  incompressibility  assumption.  This  yields  a  Poisson  equation  for  the  pressure 
fluctuations  (Chou,  1945): 
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The  solution,  given  by  Green's  theorem,  is  the  integral  of  the  right  hand  side  over  the 
entire  volume  occupied  by  the  flow  and  all  surfaces  surrounding  the  volume  (Morse 
and  Feshbach,  1953).  Thus,  the  local  pressure  fluctuations  are  determined  by  events 
taking  place  throughout  the  entire  flow  and  by  distant  boundary  conditions.  With  7^ 
as  the  current  location  and  integrating  over  all 7^,  the  pressure  fluctuations  at  7^  are 
.  .  p (!*&;  £u.du.\dV(rl) 
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surface  integrals 


(1.  10) 


Since  the  integrations  in  (1.10)  are  with  respect  to  r^,  the  pressure-rate  of 
strain  terms  can  be  formed  by  multiplying  (1.10)  by  ^u(r^  )'^Xj.  which  is  brought  into 
the  integrals,  and  then  ensemble-averaging  the  resulting  terms.  The  surface  integrals 
are  usually  neglected.  The  argument  is  that  the  influence  of  the  boundaries  is 
negligible  because  the  distance  to  the  boundaries  is  large  compared  to  the  correlation 
length  scale  for  the  pressure-rate  of  strain  terms.  This  may  be  an  inaccurate 
assumption  for  boundary  layer  flows  since  the  correlation  length  for  the  large-scale 
eddies  is  always  the  same  order  of  magnitude  as  their  distance  from  the  boundary.  The 
problem  can  also  be  recast  using  the  method  of  images  to  avoid  the  inclusion  of  the 
surface  integrals  (Launder  et  al.,  1975).  However,  it  has  not  been  clearly  demonstrated 
that  the  contribution  from  the  surface  integrals  would  be  significant.  In  this  study  the 
surface  integrals  will  be  neglected,  the  resulting  expression  is 
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where  the  primes  represent  the  expressions  evaluated  at  r2  and  the  unprimed  variables 
are  evaluated  at  The  first  integral  on  the  right  hand  side  (xn)  represents  the 
nonlinear  turbulence-turbulence  interactions.  The  remaining  volume  integral  (»trt 
represents  turbulence-mean  flow  interactions.  Unfortunately,  the  integrals  can  not  be 
solved  in  general.  However,  there  exists  a  first-order  approximation  that  is  universal; \ 


accepted  for  the  first  integral.  This  nonlinear  part  of  the  pressure-rate  of  strain  tends 
to  reduce  the  anisotropy  of  the  turbulence  (Lumley,  197S).  It  is  usually  parameterized 
with  seme  form  of  the  retum-to-isotropy  model  of  Rotta  (1951), 
nn  oc  (e  -  w^*)/t  , 

where  t  is  a  characteristic  time  scale.  Thus,  the  nonlinear  portion  of  the  pressure-rate 
of  strain  term  is  assumed  to  be  proportional  to  the  degree  of  anisotropy  of  the 
turbulence. 

The  7tr  term  represents  pressure  fluctuations  due  to  buoyancy,  mean  shear  and 
rotation,  respectively.  This  "rapid”  part  is  so  named  because  these  pressure 
fluctuations  can  respond  on  time  scales  that  are  independent  of,  and  in  principle  much 
faster  than,  the  time  scales  of  the  turbulent  velocity.  For  example,  consider  an  initially 
isotropic  turbulent  flow  that  is  subjected  to  an  intense  distortion  by  the  mean  flow. 
One  possibility  is  an  intense  velocity  shear  (Lumley,  1972).  For  some  time  after  the 
distortion,  the  turbulence  is  nearly  isotropic  so  nn  is  small  and  the  turbulence  responds 
to  only  the  mean  flow  for  some  short  time.  The  duration  of  this  period  is  determined 
by  the  intensity  of  the  distortion,  so  it  depends  only  on  the  mean  flow.  Hence,  the 
turbulence-mean  flow  terms  are  known  as  the  "rapid”  part  of  the  pressure-rate  of 
strain.  These  terms  have  been  neglected  in  most  parameterizations  of  the  pressure-rate 
of  strain  (Chapter  II). 

The  rotation  terms  in  nr  may  tend  to  cancel  the  effects  of  the  explicit  rotation 
stress  terms.  For  two-dimensional  turbulence  with  the  assumptions  that  the  pressure- 
rate  of  strain  correlation  length  scale  is  small  compared  to  variations  in  the  Coriolis 
force  and  compared  to  the  length  scale  of  the  mean  shear,  it  appears  that  these  terms 
do  cancel  the  rotation  stress  terms  (Gutteling,  1981).  However,  for  three-dimensional 
flows  this  does  not  appear  to  be  the  case  (B.  Gallagher  and  E.  Gutteling.  University  of 
Hawaii,  personal  communication). 

More  compelling  evidence  comes  from  studies  of  the  atmospheric  planetary 
boundary  layer.  Zeman  and  Tennekes  (1975)  included  a  parameterization  of  the  terms 
in  Jtr  that  involved  the  mean  shear  and  the  Coriolis  forces  in  a  model  of  the 
atmospheric  surface  (constant  stress)  layer.  If  the  turbulence  is  only  weakly 
inhomogeneous  on  scales  comparable  to  the  integral  scale,  the  mean  flow  variables  in 
7ir  can  be  taken  outside  the  integral.  The  resulting  integral  of  the  divergence  of  the 
two-point  velocity  correlations  is  a  fourth  rank  tensor.  However,  it  is  commonly 
assumed  to  be  proportional  to  the  single  point  velocity  correlations  (Launder  et  al. . 
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1975),  which  are  the  Reynolds  stresses.  The  rotation  tenns  in  the  resulting  pressure- 
rate  of  strain  parameterization  oTZeman  and  Tennekes  (1975)  are  proportional  to  the 
explicit  rotation  stress  terms.  That  is 

=  7i2cijk^juiuk  • 

Based  on  the  results  of  seven  different  shear  flow  experiments  (Table  1),  Zeman  and 
Tennekes  (1975)  concluded  that  Yj  was  0.26  ±0.04.  Results  calculated  from  the 
atmospheric  boundary  layer  experiments  were  the  two  lowest  values.  However,  this 
may  be  the  result  of  the  difficulties  in  measuring  an  environmental  flow  versus  a 
laboratory  flow,  rather  than  a  systematic  difference  between  the  atmospheric  boundary 
layer  and  other  shear  flows.  These  results  indicate  that  the  pressure-rate  of  strain 
rotation  terms  are  smaller  than  the  similar  rotation  stress  terms  that  appear  explicitly 
in  (1.1  -  1.3).  For  the  well-mixed  la>er,  the  terms  are  likely  to  be  even  smaller  than 
they  are  in  the  surface  layer  where  the  fluctuating  pressure  is  large  due  to  the  presence 
of  the  boundary.  In  this  research,  Ttr  is  parameterized  with  terms  expressed  as  small 
perturbations  of  the  explicit  terms  (see  Chapter  III  and  Garwood  et  al.,  1985a  and  b). 

TABLE  1 

NUMERICAL  VALUES  OF  Yj  FROM  ZEMAN  AND  TENNEKES  (1975) 

where  yt  is  the  coefficient  of  the  rotation  term 
in  the  rapid  part  of  their  pressure-rate  of  strain  parameterization 


Source  Yj 

Champagne  et  al.  (1970);  homogeneous  shear  flow  in  a  wind  tunnel.  0.3057 

Klebanoff  (1955);  wind-tunnel  boundary  layer  on  a  smooth  surface.  0.2928 

So  and  Mellor  (1972);  surface  layer  in  wind  tunnel,  smooth  wall.  0.2792 

Hinze  (1975);  turbulent  pipe  flow.  0.2726 

Comte-Bellot  (1965);  turbulent  flow  between  parallel  plates.  0.2644 

Cramer  (1967);  atmospheric  surface  layer.  0.2235 

Wyngaard  et  al.  (1974);  atmospheric  surface  layer,  modified  data.  0.2227 


Moeng  and  Wyngaard  (1986)  used  a  large-eddy  simulation  model  to  simulate 
the  pressure  covariance  that  results  from  the  integrals  in  (1.11)  for  a  convective 
atmospheric  surface  layer.  They  showed  that  the  pressure  covariances  due  to  rotation 
were  significantly  smaller  than  the  Jtn  terms  and  the  pressure  covariances  due  to  the 
buoyancy  fluxes.  Thus,  pressure  covariances  induced  by  rotation  in  7tr  are  probably 
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small  compared  to  the  rotation  stress  terms  and  the  other  pressure  covariances  for  the 
well-mixed  layer. 

As  indicated  above,  there  are  several  unresolved  issues  and  active  areas  of 
research  involving  the  pressure-rate  of  strain  parameterization.  However,  it  is  beyond 
the  scope  of  this  research  to  advance  the  parameterization  of  the  pressure-rate  of 
strain.  Therefore,  a  working  hypothesis  will  be  that  the  rotation  terms  in  the  pressure- 
rate  of  strain  parameterization  are  not  as  significant  as  the  explicit  rotation  stress 
terms.  Thus,  the  pressure-rate  of  strain  parameterization  used  for  this  research  will 
depend  only  weakly  on  fluctuations  due  to  mean  quantities.  There  does  exist  the 
possibility  that  any  effects  of  the  rotation  stress  terms  that  occur  in  the  numerical 
results  in  this  research  may  be  negated  by  competing  terms  in  the  pressure-rate  of 
strain.  On  the  other  hand,  if  this  research  discovers  no  significant  rotation  effects,  then 
pursuing  the  more  complex  problem  of  rotation  effects  in  the  pressure-rate  of  strain 
parameterization  is  pointless. 

C.  HYPOTHESES 

Several  working  hypotheses  will  be  tested  in  this  dissertation  to  establish  the 
effects  and  the  significance  of  rotation  stress  for  mixed  layer  dynamics.  Comparisons 
of  the  effects  of  rotation  stress  and  the  rotatior.ally  enhanced  dissipation 
parameterization  of  Garwood  (1977)  on  mixed  layer  dynamics  are  of  interest  since  both 
terms  represent  the  effects  of  rotation  on  the  TKE  budget.  However,  the  two  terms 
represent  different  physics  and  affect  the  TKE  budget  quite  differently.  Both  terms 
alter  the  rate  of  mixed  layer  deepening  and  the  equilibrium  mixed  layer  depth.  The 
rotation  stress  does  so  by  altering  the  distribution  of  TKE  between  the  horizontal  and 
vertical  components,  and  hence  altering  the  vertical  convergence  of  TKE  to  the 
entrainment  zone.  The  dissipation  enhancement  reduces  the  total  TKE  available  for 
mixing. 

The  first  hypothesis  is  that  the  functional  form  of  the  entrainment  rate  wall  be 
significantly  different  for  rotation  stress  than  for  rotationally  enhanced  dissipation. 
The  dependence  of  isotropy  of  the  TKE  budget  on  rotation  stress  should  also  be 
different  than  it  is  for  roationally  enhanced  dissipation.  Both  of  these  effects  result 
from  the  fact  that  rotation  stress  distributes  TKE  among  the  TKE  components, 
whereas  dissipation  reduces  the  TKE  isotropically.  These  hypotheses  are  tested  using 
the  nondimensional  framework  that  was  developed  in  Garwood  (1977)  and  used  to 
intercompare  mixed  layer  models  in  Garwood  (1979).  Unlike  simulations, 
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nondimensional  comparisons  are  independent  of  any  tuning  of  parameterization 
constants.  They  demonstrate  the  functional  form  of  the  models,  rather  than  the  ability 
to  fit  the  models  to  a  particular  data  set.  However,  while  nondimensional 
representations  clearly  reveal  the  funcrional  form  of  models,  they  do  not  indicate 
whether  that  functional  form  depicts  the  actual  behavior  of  the  processes  being 
modeled  or  the  interactions  among  the  processes.  To  accomplish  the  latter, 
comparisons  of  model  simulations  and  data  sets  are  required.  Therefore,  the  isotropy 
of  the  TKE  budget  and  the  deepening  rate  will  also  be  examined  using  model 
simulations. 

A  second  working  hypothesis  for  this  dissertation  is  that  diurnal  to  synoptic  scale 
shallowing  and  deepening  events  are  altered  subtantially  by  rotation  stress.  Variations 
in  wind  direction  may  substantially  change  the  equilibrium  depth  of  the  mixed  layer, 
particularly  in  late  winter  and  early  spring  when  the  mixed  layer  is  deep  and  the  surface 
wind  stress  and  surface  buoyant  damping  nearly  balance.  Both  the  frequency  of 
occurrence  and  the  vertical  scale  of  these  mixing  events  are  affected  by  rotation  stress. 
In  Garwood  ei  al.  (19S5a),  we  have  shown  that  the  equilibrium  depth  predicted  by  a 
model  that  includes  rotation  stress  is  given  by  (1.7).  This  equilibrium  depth  can  vary 
substantially  from  the  Obukhov  depth  for  east-west  surface  wind  stress  (Figure  1.3). 
This  should  be  most  apparent  at  diurnal  to  synoptic  time  scales,  w’hich  are  comparable 
to  ft  Furthermore,  it  is  argued  here  that  the  entrainment  rate  is  altered  by  the 
rotation  stress,  which  affects  deepening  events  as  well  as  shallowing  events.  While  the 
greatest  depth  which  could  be  attained  during  a  deepening  event  is  limited  by  the 
equilibrium  depth,  the  actual  depth  that  is  achieved  on  diurnal  to  synoptic  scales  is 
frequently  determined  by  the  time  scale  of  the  atmospheric  forcing.  In  those  cases,  the 
rate  of  deepening  and  the  initial  conditions,  rather  than  the  equilibrium  depth,  control 
the  vertical  scale  of  the  mixing  event. 

A  last  hypothesis  is  that  the  change  in  the  vertical  convergence  of  TKE  caused 
by  rotation  stress  may  reduce  or  eliminate  the  need  for  the  ad  hoc  rotational 
enhancement  to  the  dissipation  of  Garwood  (1977)  which  is  needed  to  prevent 
excessive  mixed  layer  deepening  for  westerly  wind  regimes.  If  rotation  stress  can 
eliminate  or  significantly  reduce  the  need  for  rotationally  enhanced  dissipation,  regions 
of  predominately  easterly  winds  would  not  have  a  one-dimensional  cyclical  steady-state 
for  mixed  layer  depth  without  a  net  downward  heat  flux,  since  rotation  stress  enhances 
mixing  for  easterly  winds.  However,  a  cyclical  steady-state  could  aiso  be  achieved  with 
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a  three-dimensional  balance  of  forces  such  as  a  mean  horizontal  advection  or 
divergence. 

The  hypotheses  will  be  tested  by  means  of  three  different  models  of  mixed  layer 
dynamics.  The  purpose  of  this  approach  is  to  distinguish  the  effects  of  various 
processes  rather  than  to  obtain  the  best  comparison  with  the  data  by  the  simulations. 
The  latter  would  probably  be  better  accomplished  with  a  combined  model  that 
incorporated  all  the  processes.  However,  because  the  interactions  between  the 
processes  are  nonlinear,  separating  the  effects  of  various  process  in  a  combined  model 
would  be  difficult  at  best. 

The  details  of  the  three  models  (ZSTAR,  HSTAR  and  RSTAR)  are  developed  in 
Chapter  III.  The  ZSTAR  model  is  the  Garwood  (1977)  model,  which  contains  the 
rotational  dissipation  enhancement  that  enters  the  equations  through  the 
nondimensional  variable  Z*.  The  HSTAR  model  is  the  model  of  Garwood  (1977) 
excluding  dissipation  enhancement  (Z*  =  0).  Thus,  HSTAR  contains  neither  the 
rotational  enhancement  to  dissipation  or  the  rotation  stress.  The  RSTAR  model  is  the 
HSTAR  model  with  the  addition  of  rotation  stress,  which  will  be  represented  in  the 
nondimensional  form  by  R*. 


II.  LITERATURE  REVIEW 

A.  BULK  OCEAN  MIXED  LAYER  MODELS 

This  brief  review  of  mixed  layer  models  will  focus  primarily  on  vertically- 
integrated  or  bulk  models.  The  purposes  are  to  place  the  ZSTAR  and  RSTAR  models 
in  historical  perspective  with  respect  to  models  of  the  same  type  and  to  review  the 
development  of  dissipation  parameterizations.  Other  types  of  models  and  closure 
schemes,  as  applied  to  the  OPBL,  will  be  mentioned  only  briefly,  and  models  of  the 
atmospheric  boundary  layer  will  not  be  reviewed  at  all.  All  the  models  being  reviewed 
use  the  TKE  equation  or  the  component  TKE  equations  and  close  the  Reynolds 
equations  at  second  order  or  lower.  Equations  for  the  shearing  stresses  are  not  used  in 
these  models.  Although  models  such  as  the  Level  3  and  4  models  discussed  in  Mellor 
and  Yamada  (1974)  do  include  equations  for  the  shearing  stresses,  this  has  not 
improved  model  performance  in  general  (Martin,  1985). 

Kraus  and  Rooth  (1961)  showed  that  an  equilibrium  could  be  established  by  a 
balance  of  surface  cooling  (due  to  longwave,  sensible  and  latent  heat  fluxes)  and 
warming  at  depth  due  to  the  absorption  of  solar  radiation.  The  equilibrium  depth  for 
this  balance  is  the  solar  radiation  compensation  depth.  Kraus  and  Turner  (1967) 
incorporated  the  mechanical  mixing  due  to  wind  stress  and  the  damping  due  to  a 
negative  surface  buoyancy  flux  (warming)  and  included  a  heat  equation.  The 
equilibrium  depth  for  this  model  is  the  Obukhov  length  scale  (L0),  which  is  the  ratio  of 
the  surface  shear  production  of  TKE  to  surface  buoyancy  flux  (1.8).  This  scale  depth 
was  first  introduced  into  ocean  modeling  by  Kitaigordsky  (1960)  using  the  atmospheric 
surface  layer  as  an  analogue.  Geisler  and  Kraus  (1969)  included  an  equation  for  the 
mean  momentum,  in  addition  to  the  TKE  and  mean  temperature  equations,  and  solved 
the  three  equation  model  for  a  slab  or  well-mixed  layer  in  which  the  mean  properties 
were  constant  (Figure  1.2). 

Miropolsky  (1970)  and  Denman  (1973)  were  the  first  to  incorporate  dissipation 
due  to  molecular  viscosity  as  a  significant  part  of  the  TKE  budget.  The  rate  of 
dissipation  (t)  is  assumed  to  depend  on  the  total  TKE  (e^  which  is  twice  the  total  TKE) 
at  the  integral  scale  (Lj)  such  that 
c  ~  e3//2  , 


29 


as  discussed  in  the  introduction.  Niiler  (1975)  added  the  loss  of  TKE  from  the  mixed 
layer  due  to  the  downward  radiation  of  internal  waves.  He  also  resolved  a  perceived 
conflict  between  turbulence  models  and  the  model  of  Pollard  et  al.  (1973).  The  later 
considered  entrainment  due  to  dynamic  instabilities  of  the  mean  current  and  used  the 
mean  kinetic  energy  rather  than  the  TKE  to  model  the  energy  balance  and  resulting 
mixing. 

Alexander  and  Kim  (1976)  and  Kim  (1976)  showed  that  an  enhancement  to  the 
dissipation  was  required  to  correctly  model  shallow  summertime  mixed  layers.  They 
hypothesized  a  constant  background  dissipation,  so  their  dissipation  parameterization 
had  the  form 

c  ■'*'  e3/2  +  Cq  , 

where  cQ  is  a  constant.  Gill  and  Turner  (1976)  showed  that  previous  mixed  layer 
models  could  not  close  the  annual  cycle  of  mixed  layer  depth  due  to  excessive 
wintertime  deepening.  Elsberry  et  al.  (1976)  used  a  dissipation  enhancement  that 
increased  exponentially  with  depth  to  improve  simulations  of  the  ocean's  response  to 
hurricanes. 

As  discussed  in  the  introduction,  Garwood  (1977)  included  a  Rossby  number 
dependent  enhancement  term  in  the  dissipation  parameterization.  He  also 
demonstrated  the  importance  of  explicitly  modeling  the  convergence  of  TKE  into  the 
entrainment  zone.  This  required  the  calculation  of  wrms!SSV'w  •  Thus,  separate  TKE 
budgets  must  be  calculated  for  the  horizontal  and  vertical  components,  whereas 
previous  bulk  mixed  layer  models  used  only  the  total  TKE  budget.  Processes  that 
transfer  energy  among  the  TKE  components  become  important  in  any  model  that 
includes  a  separation  of  the  TKE  components,  as  does  Garwood  (1977).  In  the 
introduction,  the  pressure-rate  of  strain  and  the  rotation  stress  were  shown  to  be  such 
processes. 

Garwood  et  al.  (1985a)  derived  a  bulk  model  that  included  rotation  stress  for  a 
non-entraining  mixed  layer.  We  showed  that  the  equilibrium  mixed  layer  depth  for 
that  model  is  (1.7).  In  Garwood  et  al.  (1985b),  we  computed  the  longitudinal 
variation  of  L0  and  along  the  equator  in  the  Pacific.  The  longitudinal  dependence 
of  Lc  computed  using  climatological  values  for  the  surface  heat  fluxes  and  for  the  east- 
west  surface  stress  is  shown  in  Figure  2.1.  The  vertical  bars  represent  the  variation  in 
L0  caused  by  varying  the  surface  heat  flux  by  ±  3W  m*2.  The  open  circles  are 
estimates  of  the  mixed  layer  depth  from  Lemasson  and  Piton  (1968).  The  longitudinal 
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variation  of  (Figure  2.2)  agrees  with  the  estimates  of  mixed  layer  depth 
significantly  better  than  does  LQ.  The  increased  depth  in  the  western  part  of  the  basin 
is  the  result  of  the  stronger  easterly  winds  relative  to  those  on  the  east  side  of  the 
basin. 

The  current  research  extends  the  non-entraining,  equatorial  study  of  the 
climatological,  equilibrium  mixed  layer  by  Garwood  et  al.  (1985a,  b)  to  the  case  of  the 
time-dependent,  mid  latitude  mixed  layer.  Deepening  and  shallowing  of  the  mixed 
layer  will  be  studied  for  diurnal  to  seasonal  time  scales  in  contrast  to  the  previous 
climatological  studies. 

Vertically-integrated  or  bulk  models  are  not  the  only  type  models  used  to 
simulate  the  mixed  layer.  A  variety  of  models  exist  that  use  eddy  viscosity  techniques 
to  close  the  turbulence  equations  using  only  properties  of  the  mean  field.  Since  these 
will  not  be  reviewed  here,  the  interested  reader  is  referred  to  Businger  (1982)  for  a 
general  review  of  the  equations  and  concepts  of  turbulence  modeling  in  the 
atmospheric  boundary  layer  and  a  discussion  of  eddy  viscosity  model  theory.  Another 
set  of  models  do  not  use  a  similarity  solution  to  assist  in  closing  the  equations.  Rather, 
a  set  of  turbulence  equations  is  solved  as  functions  of  z  using  standard  finite  difference 
or  finite  element  numerical  techniques.  Whereas  the  vertically-integrated  models 
require  values  only  at  the  surface  and  bottom  of  the  mixed  layer,  these  models  require 
a  set  of  dep’h-dependent  parameterizations  for  the  higher  order  moments.  A  set  of 
models  of  this  type  has  been  developed  from  the  work  of  Mellor  and  Yamada  (1974). 
This  model  was  first  applied  to  the  OPBL  by  Mellor  and  Durbin  (1975).  Other  models 
in  this  class  include  Wam-Varnas  and  Piacsek  (1979),  Kundu  (1980)  and  Klein  and 
Coantic  (1981).  A  general  review’  of  boundary  layer  modeling  using  eddy  diffusivity 
closure,  second-order  closure  and  eddy  resolving  models  is  given  by  Wyngaard  (1982). 

In  the  next  section,  studies  that  deal  with  rotating  turbulent  flows  will  be 
reviewed.  Some  background  on  pressure-rate  of  strain  parameterizations  will  be  given 
in  the  last  section. 

B.  STUDIES  OF  THE  EFFECTS  OF  ROTATION  ON  TURBULENT  FLOWS 

The  effects  of  rotation  have  been  largely  neglected  in  geophysical  boundary  layer 
modeling.  However,  the  study  of  rotation  and  of  streamline  curvature,  which  produces 
similar  effects,  has  received  considerable  attention  in  the  engineering  literature.  In  this 
section,  the  pertinent  experimental  and  theoretical  work  is  reviewed.  This  review  will 
show  experimental  evidence  for  the  existence  and  influence  of  rotation  processes. 
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particularly  rotation  stress.  The  streamline  curvature  effects  in  turbulent  boundary 
layers  that  are  equivalent  to  rotation  will  also  be  discussed.  Finally,  the  similarity  of 
engineering  models  to  the  model  proposed  in  this  dissertation  will  be  discussed. 

Two  early  works  modified  the  Richardson  number  and  mixing  length  to  include 
the  effects  of  rotation.  Rotta  (1967)  derived  a  modified  mixing  length  to  account  for 
the  effects  of  centrifugal  and  Coriolis  forces  on  TKE  production.  His  modification  was 
based  on  an  analysis  of  the  component  TKE  equations.  He  calculated  the  velocity 
distributions  for  a  turbulent  flow  between  a  fixed  cylinder  and  an  inner  concentric 
rotating  cylinder.  The  calculated  velocity  distributions  were  similar  to,  but  consistently 
less  than,  the  experimental  values  of  Taylor  (1935)  and  Wallendorf  (1935).  Bradshaw 
(1969)  extended  an  analogy,  first  proposed  by  Prantl  (1930),  between  the  effects  of 
gravity  and  the  effects  of  centripetal  or  Coriolis  acceleration.  Using  Prantl's  analogy, 
Bradshaw  constructed  Richardson  numbers  for  rotation  and  curvature  that  are 
analogous  to  the  flux  and  gradient  Richardson  numbers  for  buoyancy.  Values  for  the 
wall  shear  stress  and  the  displacement  thickness  and  momentum  thickness  of  the 
boundary  layer  calculated  from  the  curvature  Richardson  number  (centripetal 
acceleration)  compared  favorably  with  the  results  of  Schubauer  and  Klebanoff  (1951). 
However,  the  wall  shear  stress,  the  displacement  thickness  and  momentum  thickness  of 
the  boundary  layer  calculated  with  the  rotation  Richardson  number  (Coriolis 
acceleration)  were  consistently  less  than  the  experimental  values  of  Halleen  and 
Johnston  (1967).  Although  these  results  show  that  a  simple  Richardson  number  model 
is  not  completely  adequate,  such  a  model  does  yield  results  for  steady  state  cases  that 
are  consistent  with  observations.  Garwood  ei  al.  (1985b)  demonstrated  that  a  length 
scale  based  on  a  Richardson  number,  which  included  the  effects  of  rotation,  appears  to 
correctly  predict  the  longitudinal  variation  of  the  climatological  mixed  layer  depth  in 
the  equatorial  Pacific  (Figure  2.2).  This  augmented  Richardson  number  is 

h  bw  +  2fi  uw 

Ri  =  —  =  —  v 

8  Lq  \xwdU/dz 

where  is  the  rotation  length  scale  defined  in  (1.7).  This  Richardson  number 
reduces  to  the  buoyancy  Richardson  number  for  2fiyuw->0  and  to  the  rotation 
Richardson  number  for  bw-*0.  This  Richardson  number  criterion  is  derived  directly 
from  the  steady-state,  second-order  closure  model  of  this  dissertation  for  the  case  of 
zero  entrainment.  The  rotation  processes  that  transfers  energy  among  the  TKE 
components  were  explicit  terms  in  the  Garwood  et  al.  (19S5a)  model.  In  Rotta  (1967) 


and  Bradshaw  (1969).  the  rotation  processes  were  incorporated  in  the  parameterization 
of  shear  production. 

Johnston  et  al.  (1972)  reported  the  results  of  a  series  of  experiments  performed 
in  a  channel  rotating  about  an  axis  perpendicular  to  the  mean  flow  and  parallel  to  the 
side  walls  (Figure  2.3).  This  configuration  will  produce  rotation  stress  terms  analogous 
to  those  generated  in  the  OPBL.  In  the  OPBL,  east-west  winds  generate  an  east-west 
Reynolds  stress  in  the  presence  of  a  northward  rotation  vector  and  a  vertical  shear  in 
the  east-west  mean  current  (dU/dz).  In  the  experiments  of  Halleen  and  Johnston 
(1967)  and  Johnston  et  al.  (1972),  the  side  walls  generate  a  cross-stream  Reynolds 
stress  and  mean  shear  in  the  presence  of  a  vertical  rotation  vector.  A  schematic  of  the 
flow  (Figure  2.3)  shows  that  there  will  be  a  stable  and  an  unstable  side  to  the  channel, 
as  determined  by  the  direction  of  the  vorticity  generated  by  the  wall  shear  relative  to 
the  rotational  vorticity.  When  the  vorticities  are  additive,  the  flow  is  stabilized  by  the 
rotation,  and  the  wall  stress  is  reduced.  When  the  vorticities  subtract,  the  flow  is 
destabilized  by  the  rotation,  which  increases  the  wall  stress.  The  results  of  these 
experiments  showed  qualitatively  three  rotation-related  phenomena.  The  first  was  a 
change  in  the  wall-layer  streak  bursting  rate.  The  wall-layer  streak  bursting  is  due  to 
the  production  of  turbulence  at  the  wail  by  dynamic  pressure  fluctuations  in  the  inner, 
or  constant  stress,  layer.  It  is  directly  related  to  the  intensity  of  the  turbulence  in  the 
outer,  or  well-mixed  layer.  The  streak  bursting  rate  is  increased  on  the  destabilized  side 
of  the  channel  and  decreased  on  the  stabilized  side.  The  second  effect  was  the 
suppression  of  the  transition  to  turbulence  on  the  stabilized  side  of  the  channel.  The 
third  effect  was  the  formation  and  growth  of  longitudinal  roll  cells  (Taylor-Goertler 
vortices)  on  the  unstable  side.  The  results  quantitatively  showed  that  the  mixing  length 
proposed  by  Bradshaw  (1969)  provided  a  good  fit  to  the  data  for  0^Ri£0.25.  Miles 
(1961)  showed  that  laminar  flows  were  dynamically  stable  for  Ri>0.25.  Qualitative 
results  from  another  set  of  experiments  conducted  with  the  same  apparatus  (Johnston, 
1973)  indicated  that  for  small  Rossby  numbers  (rotation-dominated  flows)  the 
transition  to  turbulence  occurred  at  Reynolds  numbers  an  order  of  magnitude  larger 
than  those  for  no  rotation.  Thus,  a  stabilizing  rotation  stress  has  been  shown  to 
significantly  suppress  the  transition  to  turbulence.  However,  Johnston  (1973)  also 
found  that  the  criterion  for  transition  to  turbulence  given  by  Bradshaw  (1969)  was 
inadequate.  Thus,  the  Richardson  number  model  of  Bradshaw  does  not  account  for  all 
rotation-related  turbulence  phenomena. 
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The  literature  on  streamline  curvature  is  somewhat  more  extensive  than  for 
rotation,  probably  due  to  the  relative  ease  of  experimental  design.  Because  the  effects 
of  streamline  curvature  and  rotation  are  analogous  mathematically  (Rotta,  1967,  and 
Bradshaw,  1969),  some  experimental  and  theoretical  work  on  streamline  curvature  will 
be  included  in  this  review.  For  streamline  curvature,  centripetal  acceleration  replaces 
the  Coriolis  acceleration  of  the  rotation  case.  The  interaction  term  becomes  the 
centripetal  acceleration  multiplied  by  the  Reynolds  stress. 

lit?  _ u 

downstream  :  — — —  =  +2uw—  +  .  .  . 

2dt  R 

l5w2  _ U 

normal  to  wall:  --■■■  =  -2uw —  +  .  .  . 

2dt  R 

(see,  Rotta,  1967).  The  sign  of  the  radius  of  curvature  (R)  determines  the  direction  of 
transfer  of  TKE.  For  convex  curvature  (R>0),  the  transfer  is  from  the  component 
normal  to  the  wall  to  the  downstream  component,  which  stabilizes  the  flow. 
Conversely,  concave  curvature  (R<0)  destabilizes  the  flow  by  transferring  TKE  to  the 
component  normal  to  the  wall  (Figure  2,4). 

So  and  Mellor  (1973)  showed  that  values  for  the  critical  curvature  Richardson 
number  for  a  curved  two-dimensional  boundary'  layer  ranged  from  0.30  near  the  start 
of  the  curved  test  section  to  0.23  well  inside  the  test  section.  These  values  are  similar 
to  the  critical  buoyancy  Richardson  number  of  0.2  measured  by  Businger  et  al.  (1971) 
in  the  atmosphere  and  to  the  theoretical  value  of  0.25  (Miles,  1961).  Mellor  (1975), 
using  the  Mellor  and  Yamada  (1974)  level  2  model,  calculated  a  critical  (lux 
Richardson  number  of  0.21  and  critical  gradient  Richardson  number  of  0.22  for 
buoyancy  stabilized  flows  and  a  critical  flux  Richardson  number  of  0.19  for  the 
curvature  case.  This  is  in  good  agreement  with  experimental  results  for  the  density 
stratified  case,  but  underestimates  the  importance  of  the  curvature.  Mellor  (1975) 
hypothesizes  that  the  lack  of  turbulent  diffusion  (transport)  in  the  level  2  model  may 
account  for  the  discrepancy.  The  Garwood  (1977)  model  explicitly  includes  turbulent 
transport  through  the  vertical  convergence  of  TKE  in  the  entrainment  zone. 

Irwin  and  Smith  (1975)  used  a  modification  of  the  second-order  closure  model  of 
Launder  et  al.  (1975)  to  simulate  the  curvature  stress  effects  found  in  the  experiments 
of  Gillis  et  al.  (1965),  Schubauer  (1951)  and  Meroney  and  Bradshaw  (1975).  In  that 
model,  the  curvature  stress  was  parameterized  as  an  extra  production  term  (the  "extra 
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rate  of  strain"  concept  of  Bradshaw.  1973).  Hunt  and  Joubert  (1979)  used  a  second- 
order  closure  model  that  explicitly  accounted  for  the  "conservative  reorientation"  of.  or 
the  transfer  of  TKE  between,  the  TKE  components  by  the  curvature  stress  to  explain 
the  results  of  their  experiments  with  "small"  streamline  curvature  (large  R).  Here, 
small  streamline  curvature  implies  that  the  curvature  does  not  significantly  affect  the 
mean  flow.  Gibson  et  al.  (1981)  simulated  the  response  of  a  flow  to  large  curvature. 
They  used  a  second-order  closure  model  with  the  unsteady  terms  parameterized  in 
terms  of  transport  by  the  mean  flow,  and  the  pressure-rate  of  strain  parameterization 
included  approximations  for  "rapid"  distortion  effects  (rapid  distortion  will  be  discussed 
in  the  next  section).  The  model  was  able  to  simulate  the  results  of  Gillis  and  Johnston 
(1979)  for  the  curved  region,  but  not  for  the  region  of  recovery  to  a  fiat  boundary 
layer. 

Muck  et  al.  (1985)  studied  the  response  of  a  fully  developed  turbulent  flow  to  a 
sudden  convex  curvature.  In  a  companion  paper,  Hoffman  et  al.  (1985)  studied  the 
response  of  a  fully  developed  turbulent  flow  to  a  sudden  concave  curvature.  They 
concluded  that  the  turbulent  responses  to  convex  and  concave  curvature  were 
fundamentally  different  despite  the  similarity  of  the  dimensional  analysis  for  the  two 
cases.  Their  argument  is  based  on  the  difference  in  the  structure  of  the  turbulence. 
For  concave  curvature,  the  centripetal  acceleration  destabilizes  the  flow  and  Taylor- 
Goertler  vortices  are  formed,  as  in  the  rotation  experiments  of  Johnston  et  al.  (1972). 
For  convex  curvature,  the  centripetal  acceleration  stabilizes  the  flow  and  the  turbulence 
is  damped,  but  its  structure  is  virtually  unaffected.  This  is  analogous  to  the  problem  of 
stable  versus  unstable  buoyancy  regimes.  Businger  et  al.  (1971)  empirically  determined 
a  stability  correction  to  the  mixing  length  in  the  atmospheric  surface  layer.  They  found 
that  the  form  of  the  stability  correction  as  a  function  of  nondimensional  height  was 
different  for  stable  versus  unstable  surface  boundary  layers.  The  significant  variations 
in  the  structure  of  the  turbulence  between  stable  and  unstable  regimes  and  the  effects 
of  these  structural  changes  deserves  future  study.  For  this  research,  the  flow’  structure 
will  be  assumed  to  be  similar  in  both  the  unstable  and  stable  regimes,  whether  due  to 
rotation  or  buoyancy.  This  assumes  that  the  flow  is  dominated  by  the  shear,  so  that 
the  structural  differences  in  the  turbulence  associated  with  the  unstable  and  stable 
regimes  are  not  significant. 

The  rotation  stress,  or  the  equivalent  streamline  curvature  effect,  does  not  affect 
every  rotating,  or  curved,  turbulent  flow.  The  rotation  stress  can  be  written  in  a 
generalized  form  as 
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where  Cjjk  *=  0  if  any  two  indices  are  repeated  (Garwood  ei  al.,  1987).  Fernando  (19S7) 
mentioned  several  experiments  for  which  the  rotation  stress  is  zero  because  the  rotation 
axis  is  not  orthogonal  to  the  Reynolds  stress  (i.e.,  two  indices  repeated  in  the  above 
expression).  There  must  be  a  component  of  the  Reynolds  stress  in  the  plane 
perpendicular  to  the  axis  of  rotation  for  a  nonzero  rotation  stress.  Conversely,  rotation 
stress  is  not  the  only  mechanism  by  which  rotation  affects  turbulent  flows.  Other 
mechanisms  such  as  inertial  waves  and  centripetal  acceleration  can  also  play  significant 
roles. 

Ibbetson  and  Tritton  (1975)  found  that  the  decay  rate  of  grid-stirred  turbulence 
in  a  rotating  tank  increased  with  increasing  rotation  rate.  Also,  the  ratio  of  the  parallel 
to  the  perpendicular  turbulence  length  scales  relative  to  the  rotation  axis  increased  with 
increasing  rotation  rate.  However,  the  ratio  of  turbulent  intensities  in  the  two 
directions  remained  essentially  constant.  In  this  experiment,  the  turbulence  was 
generated  by  oscillating  a  grid  parallel  to  the  rotation  axis.  No  rotation  stress  was 
produced  since  there  was  no  stress  perpendicular  to  the  rotation  axis.  In  this  case,  the 
increased  decay  rate  is  probably  due  to  the  generation  of  inertial  waves  and  the  change 
in  length  scales  is  probably  the  result  of  centripetal  accelerations. 

Hopfinger  ei  al.  (19S2)  used  a  grid  oscillating  parallel  to  the  axis  of  rotation  in  a 
rotating  tank  to  study  the  transition  from  three-dimensional  turbulence  near  the  grid  to 
two-dimensional  turbulence  away  from  the  grid.  Again,  the  rotation  stress  is  zero. 
Away  from  the  grid,  the  scale  of  the  turbulence  increases  and  the  turbulent  intensity 
decreases,  so  that  the  rotation  effect  becomes  increasingly  important.  For  the  Rossby 
number  ~ 0.2,  the  flow  breaks  down  into  a  two-dimensional  field  of  vortices  aligned 
with  the  rotation  axis  and  remains  in  this  configuration  with  increasing  distance  from 
the  grid.  This  pattern  is  consistent  with  the  results  of  Ibbetson  and  Tritton  (1975) 
regarding  the  change  in  length  scale  with  time.  A  similar  set  of  experiments  by 
Dickinson  and  Long  (1983)  gave  similar  results.  Maxworthv  (1986)  used  the  apparatus 
of  Hopfinger  et  al.  (1982)  to  show  that  the  entrainment  rate  across  a  density  interface 
decreases  with  decreasing  Rossby  number  (increasing  rotation  rate).  Based  on  the 
entrainment  model  of  Linden  (1973),  which  did  not  consider  the  effects  of  rotation. 
Maxworthv  tentatively  concluded  that  the  change  in  entrainment  rate  with  rotation 
was  due  to  the  loss  of  TKE  to  inertial  oscillations  induced  in  the  mixed  layer.  Again, 
no  rotation  stress  was  generated  with  the  configuration  used  in  this  apparatus. 
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In  summary,  the  experimental  evidence  shows  significant  effects  due  to  rotation 
stress  and  to  curvature  stress.  The  attempts  to  simulate  these  processes  numerically 
indicate  that  a  second-order  closure  model  that  explicitly  accounts  for  the  conversion 
of  TKE  among  components  is  necessary  to  model  the  phenomenon  of  rotation  stress. 
However,  the  appropriate  augmented  Richardson  number  can  give  approximate  results 
for  the  steady-state  case. 

C.  PRESSURE-RATE  OF  STRAIN  PARAMETERIZATIONS 

As  discussed  in  the  introduction,  it  is  possible  that  rotation  effects  in  the 
pressure-rate  of  strain  parameterization  may  counteract  the  explicit  rotation  stress. 
However,  the  exact  effects  of  the  various  processes  that  contribute  to  the  pressure 
fluctuations  are  not  well  documented  by  experiment  nor  are  they  well  understood 
theoretically.  Thus,  parameterizing  the  pressure-rate  of  strain  is  difficult  and 
controversial.  This  review  will  focus  on  some  aspects  of  the  parameterization  of  the 
the  rapid  part  of  the  pressure-rate  of  strain  (Jtr). 

Parameterizations  for  Jtr  can  largely  be  divided  into  two  categories:  the  first 
assumes  a  specific  structure  for  the  turbulence;  the  second  assumes  that  the  integral 
can  be  expressed  as  a  series  expansion  about  an  isotropic  TKE  distribution  (Wyngaard, 
1982).  The  difficulty  with  the  former  is  that  the  parameterization  can  only  be  valid  for 
geometrically  similar  flows,  and  the  geometries  chosen  must  be  fairly  simple  or  the 
integral  is  still  not  solvable.  The  difficulty  with  the  second  is  that  important  processes 
that  these  terms  are  supposed  to  represent  (e.g.,  rapid  distortions)  are  not  a 
perturbation  on  an  isotropic  state.  Both  methods  neglect  the  surface  integrals  in  the 
general  expression  for  p.  Thus,  they  ignore  wall  effects  that  are  fundamental  to 
boundary  layer  flows.  Some  investigators  (Gibson  and  Launder,  1978)  have  included 
boundary  effects  in  a  parameterization  of  the  pressure-rate  of  strain,  but  the  results  do 
not  seem  to  constitute  a  measurable  improvement  (Wyngaard,  1980). 

Donaldson  (1972)  discussed  a  set  of  theoretical  criteria  for  the  parameterization 
of  second  and  higher  order  terms  in  the  Reynolds  averaged  equations.  To  satisfy  these 
constraints,  the  parameterization  must  preserve  the  tensor  character  and  symmetry 
properties  of  the  original  terms,  be  dimensionally  correct,  be  invariant  under  Galilean 
transformations  and  must  satisfy  the  conservation  laws  that  were  satisfied  by  the 
original  terms.  The  fact  that  no  parameterization  of  the  pressure-rate  of  st  lin  has  yet 
been  proposed  that  can  meet  all  the  theoretical  constraints,  and  not  be  so  cumbersome 
as  to  be  unusable,  points  out  the  difficulties  of  parameterizing  the  pressure-rate  of 


strain.  Schumann  (1977)  further  noted  the  importance  of  realizability  for  turbulence 
parameterizations.  Realizability  requires  that  the  parameterization  to  produce  only 
values  that  can  be  realized  in  actual  flows.  That  is.  the  model  should  not  be  able  to 
generate  negative  TKE,  negative  TKE  components,  negative  dissipation  or  covariances 
greater  than  one  (Lumley,  1978).  Lumley  (1979)  attempted  to  construct  a  model  that 
met  these  requirements,  but  the  resulting  expressions  were  too  cumbersome  to  be  cf 
practical  use  (Wyngaard,  1982). 

The  parameterization  of  the  pressure-rate  of  strain,  particularly  the  rapid  part  of 
the  term,  has  remained  a  major  problem  in  turbulence  closure  (see  reviews  by 
Reynolds,  1976,  Lumley,  1979  and  Zeman,  1981).  Recent  work  has  attempted  to 
obtain  more  universal  and  fundamental  expressions  for  the  pressure-rate  of  strain  by 
limiting  the  number  and  type  of  assumptions  used  in  the  parameterization.  As 
discussed  in  Chapter  I,  the  integral  of  the  two-point  correlations  of  the  velocity 
derivatives  that  results  from  the  expression  for  the  rapid  part  of  pressure-rate  of  strain 
(1.11)  is  a  fourth  rank  tensor.  Gallagher  ei  al.  (1981)  started  with  a  general  expression 
for  this  fourth  rank  tensor,  which  is  the  sum  of  81  terms  each  with  a  scalar  coefficient, 
and  by  rigorous  arguments  reduced  the  number  of  unknown  coefficients  from  81  to  6. 
The  six  remaining  coefficients  required  assumptions  regarding  the  nature  of  the 
particular  flow.  Morris  (1984)  obtained  a  general  form  for  the  pressure-rate  of  strain 
parameterization  that  was  similar  to  that  of  Launders  ei  al.  (1975)  without  requiring 
that  the  turbulence  be  homogeneous.  This  research  is  not  intended  to  obtain  a  general 
solution  to  the  pressure-rate  of  strain  closure  problem.  However,  the  possible  effects 
of  the  pressure-rate  of  strain  on  the  TKE  budget  and  on  the  rotation  stress  terms  must 
be  recognized. 

;  In  summary,  this  review  has  shown  that  most  models  of  geophysical  mixed  layers 

ft  have  been  constructed  in  such  a  way  that  rotation  processes  are  neglected.  Models 

that  do  alio"''  for  rotation  effects  have  been  developed  for  other  flow  applications  and 
have  improved  predictive  capability  for  rotating  boundary  layers  with  respect  to  models 
that  neglect  rotation.  Also,  there  is  a  substantial  body  of  evidence  from  the 
engineering  literature  that  rotation  plays  an  important  role  in  turbulent  boundary 
layers  that  are  dynamically  similar  to  the  OPBL.  Thus,  there  is  reason  to  expect  that 
including  rotation  stress  in  a  model  of  the  OPBL  should  improve  the  predictive 
capability  of  the  model.  Furthermore,  the  investigation  of  a  more  complete  and 
sophisticated  picture  of  the  OPBL  is  worthwhile  on  its  own  merits. 
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III.  THEORY 


A.  MODEL  DEVELOPMENT 

The  model  developed  in  this  research  assumes  vertical  homogeneity  for  the  mean, 
dynamical  variables  for  -h  ^  z  ^  -6W.  The  well-mixed  profile,  discussed  in  Chapter  I 
(Figure  1.2),  is  consistent  with  the  average  state  of  the  upper  ocean.  The  value  of  h  is 
assumed  to  correspond  closely  to  the  mixed  layer  depth,  and  "mixed  layer"  here  is 
synonymous  with  turbulent  boundary  layer.  However,  h  is  formally  the  depth  of  the 
maximum  vertical  penetration  of  TKE  and  turbulent  fluxes  generated  in  response  to 
the  wind  stress, 

r=  'A  +  tA  *  "  Puw(°)£x  -  pw(0)£y  ,  (3.1) 

and  the  surface  buoyancy  flux,  bw(0). 

The  model  described  in  this  chapter  is  for  a  time-dependent  mixed  layer  that  may 
deepen  or  retreat  in  response  to  the  surface  fluxes.  The  mixed  layer  deepens  in 
response  to  increased  TKE  production  by  wind  stress  or  surface  cooling  (positive 
surface  buoyancy  flux)  when  water  is  entrained  into  the  fully  turbulent  OPBL  from  the 
pycnocline  below.  When  the  mixed  layer  retreats,  a  new  shallower  mixed  layer  forms. 
The  portion  of  the  water  column  that  had  been  part  of  the  mixed  layer  now  contains 
water  of  almost  homogeneous  properties,  but  the  flow  is  no  longer  actively  turbulent. 
The  following  discussion  parallels  section  2  of  Garwood  et  al.  (1985a),  who  discussed 
the  equilibrium  or  time-independent  simplification  of  the  more  general  theory' 
developed  in  this  research. 

Integrating  (1.1),  (1.2)  and  (1.3)  through  the  mixed  layer  yields: 

0  =  Gx  +  nxx  +  Qzh<uv>  +  Qvh<uw>  -  D/3,  (3.2) 

0  =  Gy  +  nyy  -  Hzh<mT>  -  D/3,  (3.3) 

0  =  <Tow>  +  nzz  -  nyh<uw>  -  D/3,  (3.4) 

where 


d  <u2> 
dt  2 

d  < v^> 

at  2 

a  <vr> 
at  2 


dz  = 


dz  = 


dz  = 
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<  >  =  dz 

-H-f 

These  equations  are  derived  assuming  local  isotropy  at  the  Kolmogorov  length  scale, 
linear  flux  profiles,  and  one  dimensionality  for  the  TKE  budget  (a  good  approximation 
even  if  the  density  and  momentum  fields  are  not  horizontally  homogeneous,  see 
Chapter  I).  We  assume  a  steady-state  TKE  budget  and  set  the  left  hand  side  of 
(3.2)-(3.4)  to  zero  since,  as  discussed  previously,  the  unsteady  terms  are  important  only 
on  very  short  time  scales  (less  than  100  seconds).  The  integral  through  the  OPBL  of 
the  shear  production  plus  the  production  due  to  breaking  surface  waves  minus  the  loss 
due  to  the  downward  radiation  of  internal  waves  (terms  G.l  -  G.3  and  P.4  in  1.1  -  1.3) 
is  G.  The  vertically-integrated  pressure-rate  of  strain  is  II,  (11.1  -  11.3  in  1.1  -  1.3).  The 
viscous  dissipation  rate  (£)  in  (1.1  -  1.3)  integrated  through  the  mixed  layer  is  D.  The 
terms  (R.l  -  R.4  in  1.1  -  1.3),  namely 

are  the  rotation  stresses  (the  interaction  of  rotation  of  the  earth  with  the  Reynolds 
shear  stresses).  The  sign  of  <uw>  is  such  that  rotation  stress  increases  uz  in  (3.2)  if 
<Tiw>  is  positive.  The  planetary  rotation  vector  is 

si  =  +  n  & 

y  y  z  z 

where  |llj  =  7.29  x  lO^sec'1.  For  the  mean  horizontal  momentum  to  be  independent 
of  depth,  the  vertical  turbulent  momentum  fluxes  (uw  and  vw)  in  the  mixed  layer  must 
be  linearly  dependent  on  z  or  constant.  If  the  depth  dependence  is  assumed  to  be 
linear,  the  vertically  integrated  momentum  fluxes  (the  Reynolds  shear  stresses)  are 

_  h  _  _ 

<uw>  =  -(uw(-h)  -  uw(0))  , 

2 

_  h  _  _ 

<vw>  *  -(vw(-h)  -  vw(  0) )  , 

2 

where  uw(0)  and  vw(0)  are  defined  in  (3.1).  The  fluxes  at  the  base  of  the  mixed  layer 
are  given  by  the  "jump  conditions" 

,  and 


where 
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uv7(  -h)  = 
W(  -h)  = 


AU  =  <U>  -  U(  -h-5 ) 


and 

AV  =  <V>  -  V( -h-6) 

are  the  mean  velocity  change  through  the  entrainment  zone  (see  Niiler,  1975). 

The  focus  of  this  study  is  on  the  exchange  of  energy  between  u-  and  w2  by  the 
action  of  the  east-west  wind  stress  (rx  =  -uw(0))  in  the  presence  of  the  northward 
component  of  planetary  rotation  (fiy).  This  term  will  reduce  vertical  mixing  if  the 
winds  are  westerly  (t  >0),  and  it  will  increase  vertical  mixing  if  the  winds  arc  easterly 
(tx<0). 

The  total  shear  production  plus  the  production  of  TKE  by  breaking  surface 
waves  minus  the  loss  at  z  =  -h  due  to  the  downward  radiation  of  internal  waves 
integrated  over  the  mixed  layer  is 


The  parameterization  for  G  is 

|tl3/2  [< AU)2+(  AV)2]3h 

G  =  m ,  “*  +  • 1  — 1,1  ■  , 

3  p  2  3t 

where  m3  and  the  subsequent  m,,  m2,  m4  and  m5  are  dimensionless  tuning  coefficients. 
These  coefficients  are  assumed  to  be  reasonably  constant  for  the  upper-ocean 
boundary  layer  under  a  wide  variety  of  atmospheric  forcing  conditions.  The 
production  of  TKE  by  the  shear  of  the  mean  current  will  be  neglected  in  this  study. 
This  may  not  be  a  good  assumption  when  h  is  small  or  in  regions  of  strong  surface 
currents.  Martin  (1985)  speculates  that  mixed  layer  simulations  of  the  shallow 
summertime  mixed  layer  might  be  improved  by  the  inclusion  of  TKE  production  by  the 
mean  shear  in  the  entrainment  zone.  However,  deSzoeke  and  Rhines  (1975)  and 
Garwood  and  Yun  (1979)  demonstrate  that  wind  shear  production  of  TKE  is  the  more 
important  source  of  energy  for  deepening  the  mixed  layer  in  general.  Therefore.  TKE 
production  by  the  shear  of  the  mean  current  will  be  neglected  for  this  study,  which  is 
primarily  concerned  with  late  winter  and  early  spring  cases  for  which  the  mixed  layer  is 
deep.  The  ratio  of  the  integral  of  the  rotation  stress  associated  with  the  mean  shear 
through  the  mixed  layer  to  that  of  the  production  of  TKE  by  the  mean  shear, 
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n  hAU— 
y  dt 


20  h 


£(  AU)~  +  (AV)2]dh  AU  ' 

2  dt 

is  only  O(10'2)  even  for  a  1  m  s"1  shear  in  the  mean  current  and  a  100  m  deep  mixed 
layer.  Therefore,  the  rotation  stress  associated  with  the  mean  shear  will  also  be 
neglected. 

For  the  mean  buoyancy  to  be  independent  of  depth,  the  buoyancy  flux  in  the 
mixed  layer  must  be  linearly  dependent  on  z  or  constant.  The  vertically  integrated 
buoyancy  flux  is 

, _  h _ 

<Bw>  =  -(bw(-h)  -  u*b* )  . 

The  effective  buoyancy  flux  through  the  surface  (u*b*)  is  the  net  downward  surface 
buoyancy  flux  due  to  the  sensible  and  latent  heat  fluxes  and  longwave  radiation  plus 
the  buoyancy  flux  due  to  the  absorption  of  solar  radiation  (Qs)  as  a  function  of  z: 


u*b*  =  -bw( 0 ) 


—  JlQs  - 

Pcp  -K-r 


2 

h+6 


i  dz 


The  absorption  of  solar  radiation  with  depth  significantly  affects  mixed  layer  dynamics, 
particularly  for  shallow  summertime  mixed  layers  (Martin,  1985;  Simpson  and  Dickey, 
1981).  The  flux  of  buoyancy  through  the  base  of  the  mixed  layer  is  due  to  the 
entrainment  of  the  more  dense,  underlying  water.  This  entrainment  buoyancy  flux  is 


bw(-h)  =  - 


ABdh 

dt 


where  AB  =  <B>  -  B(-h-8)  is  the  mean  buoyancy  change  through  the  entrainment 
zone. 


The  total  dissipation, 

f*£dz  =  D  , 

-V.-f 

is  a  function  of  the  vertical  average  of  the  turbulent  kinetic  energy, 


D  =  m^u-yV2 

where  "oil  =*  u2  +  v^V  (i.e.,  the  normal  summation  convention  on  repeated  indices 
is  implied).  Although  the  large-scale  turbulence  is  anisotropic,  the  production  and 
dissipation  scales  are  assumed  to  be  sufficiently  well  separated  (i.e.,  large  Reynolds 
number)  so  that  the  turbulence  is  isotropic  at  the  dissipation  (or  Kolmogorov)  scale. 
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This  implies  that  the  dissipation  should  be  partitioned  equally  among  the  normal 
Reynolds  stress  components.  For  the  ZSTAR  model,  the  dissipation  parameterization 
is 


D„  =  m,(  <^ou>1//2  +  — ^-fh)  <u.u;>  .  (3.5) 

«  1,1  mj  1 

where  the  additional  term  is  the  rotational  dissipation  enhancement.  By  definition,  the 
rotation  stress  is  neglected  in  the  ZSTAR  model. 

The  extent  of  vertical  mixing  depends  upon  the  rate  of  conversion  of  horizontal 
turbulent  kinetic  energy  (  <  u2  +  v2  >  )  to  vertical  turbulent  kinetic  energy  (  <  w2  >  ).  As 
discussed  in  Chapter  I,  there  are  two  processes  that  may  accomplish  this.  The  first 
process  is  the  pressure-rate  of  strain  given  by 

n  =  iiZE 

««  p  dxa 

with  no  summation  implied  on  repeated  indices.  The  second  process  is  the  rotation 
stress, 

-VWv'p  • 

The  pressure-strain  rate  correlation  (naa)  can  generally  be  divided  into  two 
parts:  the  nonlinear  self-interaction  of  the  turbulent  field  (n*\a)  and  the  "rapid"  terms 
resulting  from  the  interaction  of  the  turbulence  and  the  mean  field  (nRaa),  as 
discussed  in  Chapter  II.  The  n^aa  terms  cause  a  reduction  in  the  anisotropy  of  the 
integral  scale  turbulence  and  are  conventionally  parameterized  using  the  return-to- 
isotropy  model  of  Rotta  (1951).  A  version  of  this  parameterization,  which  is  consistent 
with  our  assumptions,  is  _ 

nNaa  *  m2<u^;>l/2(  <uju*>  -  3<ua2>). 

The  "rapid"  terms  were  discussed  at  length  in  Chapters  I  and  II.  As  discussed  in 
Chapter  II,  this  research  is  not  intended  to  obtain  a  general  solution  to  the  pressure- 
rate  of  strain  closure  problem.  Rather,  a  simple  approximation  for  IIRaa  will  be  used 
and  this  study  will  focus  on  the  impact  of  the  explicit  rotation  stress  term.  The 
approximation  that  will  be  used  for  nRoa  has  been  used  effectively  to  model  a  variety 
of  turbulent  flows  that  do  not  depend  strongly  on  "rapid"  distortion  effects. 

Expressions  for  IIRtta  can  be  defined  in  terms  of  the  vertically  integrated 
variables  discussed  above.  A  set  of  such  parameterizations  is: 


nR  =  - 

XX 


alGx 


-  a2fizh<uv>  +  a3fly<uw> 


and 


nRyy  =  -  a4Gy  -  a2fizh<uv>  , 

nRzz  =  -  a5<bw>  -  a3Hy<uw>  . 

These  are  consistent  with  the  work  of  Zeman  and  Tennekes  (1975),  which  was 
discussed  in  Chapter  I.  Thus,  the  rapid  part  of  the  pressure  terms  can  be  incorporated 
into  the  other  turbulence-mean  flow  terms.  The  literature  review  in  Chapter  II 
indicated  the  uncertain  state  of  knowledge  of  the  magnitude  and  effects  of  these  terms 
in  general.  In  particular,  there  is  a  lack  of  experimental  evidence  concerning  the 
magnitudes  and  effects  of  the  nRaa  on  geophysical  turbulent  boundary  layers.  In 
view  of  this  situation,  we  will  assume  that  the  a^  <  <  1  for  all  i.  This  is  consistent  with 
Zeman  and  Tennekes  (1975),  who  estimated  that  these  rotation  terms  were  small 
(a,  =  0.26 ±0.04)  for  the  atmospheric  surface  layer  (See  Chapter  I).  This  is  also 
consistent  with  the  results  of  Zeman  and  Lumley  (1976)  and  Moeng  and  Wyngaard 
(1986).  Also,  Wyngaard  et  al.  (1974)  evaluated  the  effect  of  the  explicit  Coriolis  terms 
on  the  structure  of  the  Reynolds  stress  profiles  for  a  buoyancy  driven  atmospheric 
planetary  boundary  layer  and  concluded  that  there  was  no  significant  change  in  the 
structure  of  the  profiles.  The  present  study  concerns  the  stress-driven  ocean  planetary 
boundary  layer.  However,  the  assumption  that  the  Reynolds  stress  profiles  are  linear 
in  z  is  consistent  with  these  results  of  Wyngaard  et  al.  (1974). 

B.  NONDIMENSIONAL  EQUATIONS 

The  TKE  equations  can  be  written  in  nondimensional  form  to  intcrcompare 
various  models  and  model  terms.  Normalizing  the  terms  in  (3.2  -  3.4)  with  the  surface 
shear  production  yields  the  nondimensional  variables  given  in  Table  2.  The  RSTAR 
model  developed  in  this  research  includes  rotation  stress  and  neglects  dissipation 
enhancement.  The  nondimensional  equations  are  shown  in  Table  3.  The  dependent 
variables  are:  P*.  the  nondimensional  entrainment  rate;  E*,  the  nondimensional  total 
TKE;  and  W*2,  the  vertical  component  of  the  TKE.  The  independent  variables  are  the 
nondimensional  mixed  layer  stability  (H*)  and  the  nondimensional  surface  rotation 
stress  (R*).  Equations  (3.7a)  -  (3.7c)  form  a  nonlinear,  coupled  set  of  equations  for 
E*,  P*  and  W*2  as  functions  of  H*  and  R*. 

Two  other  models  will  be  used  for  comparison.  They  are  the  ZSTAR  model 
(Table  4)  from  Garwood  (1977)  and  the  HSTAR  model  (Table  5).  The  ZSTAR  model 
contains  a  dissipation  enhancement  term  (Z*),  which  depends  on  rotation  but  does  not 
contain  rotation  stress  (R*).  Garwood  (1977)  demonstrated  that  the  dependence  of  P* 


TABLE  2 

THE  NONDI MENSIONA L  VARIABLES 


hAB  dh 
2m3u*  dt 


entrainment  buoyancy  flux 
wind  shear  production 


( 3. 6a) 


total  TKE _ 

surface  momentum  flux 

vertical  TKE 
surface  momentum  flux 


( 3. 6b) 


(3.  6c ) 


mixed  layer  depth 
rotation  length  scale 


( 3. 6d) 


H* 


-bw(  0)h 
2m3u*^ 


surface  buoyancy  flux 
wind  shear  production 


( 3. 6e ) 


R* 


nvhtx/P 

J  '  f  / 
2m3u* 


surface  rotation  stress 

-  ;  (  3. 6f ) 

wind  shear  production 


where  Pj 


p2  =  mj/m,,  p3  =  m5/m, 
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TABLE  3 

THE  RSTAR  MODEL:  P*  =  P*(R*.  H*) 


2E*3/2 

3 

p,E*2(E*  -  3W*2)  +  R*  , 

( 3. 7a) 

e*3/2 

0  =  -H*  -  P*  -  —  + 

3 

p2E*2(E*  -  3W*2)  -  R*  , 

( 3. 7b) 

p,W*E* 

0  =  -P*  +  -i  --  . 

2 

TABLE  4 

(  3.  7c) 

THE  ZSTAR  MODEL:  P*  =  P*(Z*,  H*) 


2E*3/2 

0  =  1  -  (E*  +  Z*) 

-  p2E*2(E*  -3W*2), 

(  3.  8a) 

E  3/2 

0  =  -H*  -  P*  -  -i-  (E*  +  Z*) 

o 

+  p2E*2(E*  -  3W*2), 

(3.8b) 

p.  W*E* 

0  =  -P*  +  H  . 

2 

(3.  8c) 

either  a  negative  net  buoyancy  flux  or  advection.  However,  the  physical  basis  for  this 
enhanced  dissipation  is  unclear,  and  Z*  may  represent  physical  processes  that  are  net 
actually  due  to  the  dissipation.  These  effects  may  include  upwelling  due  to  large-scale 
horizontal  divergence,  loss  of  TKE  in  the  mixed  layer  due  to  internal  waves  (Muller  et 
al.,  I9S4),  inertial  oscillations  and  rotation  stress. 

The  variables  in  (3.Sa)  -  (3.8c)  from  Table  4  are  the  same  as  in  (3.7a)  -  (3.7c), 
except  with  the  addition  of  rotationally  enhanced  dissipation  Z*,  defined  by  (3.6d), 
which  is  the  inverse  of  a  boundary  layer  Rossbv  number  and  the  elimination  of  the 
rotation  stress  R*.  The  dissipation  in  the  ZSTAR  model  is 

D  =  E*(E*1/2  +  Z*)  , 

which  is  the  nondimensional  equivalent  of  the  enhanced  dissipation  (3.5).  In  both  the 
ZSTAR  model  and  the  RSTAR  model,  the  dissipation  is  distributed  equally  among  the 
three  TKE  components.  The  term  E*Z*  imposes  an  additional  length  scale  (f'u*)  that 
affects  h.  The  equilibrium  mixed  layer  depth  is  now  a  function  of  a  boundary  layer 
Richardson  number  and  a  boundary  layer  Rossby  number.  Expressed  in  terms  of 
length  scales,  this  yields  h  =  G(L0,Vu*/f).  This  is  similar  to  the  equilibrium  depth 
scale  for  rotation  stress  (1.7)  that  was  introduced  in  Chapter  1. 

For  Z*  =  R*  =  0,  the  two  models  RSTAR  and  ZSTAR  are  identical.  The 
resulting  model  (HSTAR)  depends  on  H*  (i.e.  surface  buoyancy  flux)  exactly  as  the 
RSTAR  and  ZSTAR  models  do  (3.9b  from  Table  5).  However,  the  HSTAR  model 
does  not  include  rotation  stress  or  a  dissipation  enhancement.  Thus,  comparisons  of 
the  RSTAR  and  ZSTAR  models  with  the  HSTAR  model  allow  the  effects  of  the 
rotation  stress  and  the  dissipation  enhancement  respectively  to  be  isolated. 

The  terms  R*  and  Z*  are  the  inverses  of  two  boundary  layer  Rossby  numbers, 
which  represent  the  ratio  of  the  shear  production  of  TKE  to  the  vertical  and  horizontal 
rotation  stresses  respectively.  Although  the  mathematical  forms  of  R*  and  Z*  are 
similar,  the  processes  impact  the  TKE  budget  in  fundamentally  different  ways.  The 
energy  redistribution  (R*)  represents  a  known  physical  process  that  results  from  vortex 
stretching  by  the  interaction  of  the  Reynolds  stresses  and  planetary  rotation.  It  is  a 
direct  integration  of  the  original  terms  that  involves  no  higher  order  moments  and 
requires  no  parameterization.  Conversely,  Z*  is  part  of  a  parameterization  of  the 
dissipation  of  turbulence  by  molecular  viscosity.  As  discussed  in  Chapter  II.  the 
parameterization  of  dissipation  in  mixed  layer  models  has  become  increasingly 
sophisticated.  However,  it  has  usually  been  assumed  that  TKE  is  dissipated 


isotropically.  As  will  be  shown  in  the  following  sections,  this  assumption  may  be 
incorrect  for  dissipation  in  geophysical  planetary  boundary  layers. 

Since  the  independent  variables  R*  and  Z*  enter  the  equations  with  different 
powers  of  E*,  the  effects  on  the  mixed  layer  ought  to  be  dissimilar.  In  RSTAR 
compared  to  ZSTAR,  the  equations  have  been  simplified  by  eliminating  one  of  the 
highest  order  terms  in  E*,  namely  E*3  2Z*.  The  number  of  constants  also  is  reduced 
by  one,  since  p3  in  the  definition  of  Z*  is  eliminated.  At  the  same  time,  the  number  of 
physical  degrees  of  freedom  is  increased  since  wind  direction  now  affects  the  turbulent 
boundary  layer  dynamics.  Mathematically,  this  means  that  R*  less  than  zero  is 
allowed  whereas  Z*  less  than  zero  is  not.  This  will  be  discussed  in  more  detail  in  the 
next  section. 

C.  IMPLICATIONS  FOR  THE  ENERGY  BUDGET 

Most  integral  models  of  the  OPBL  are  based  on  the  total  TKE  equation.  The 
models  of  Kraus  and  Turner  (1967;,  Niiler  (1975),  Elsberry  and  Camp  (1978)  and  Kim 
(1976)  are  examples.  For  these  models,  the  rate  of  mixed  layer  deepening  is  determined 
completely  by  the  boundary  conditions  on  the  mixed  layer.  However,  Garwood  (1977) 
demonstrated  the  importance  of  the  convergence  of  TKE  at  z  =  -h  for  mixed  layer 
deepening.  Thus,  both  the  horizontal  (u2  +  v2)  and  the  vertical  (w~)  TKE  components 
must  be  known.  As  a  result,  the  transfer  of  TKE  between  the  horizontal  and  vertical 
components  becomes  important,  and  the  sources  and  sinks  of  TKE  must  be  considered 
in  relation  to  the  TKE  component(s)  they  directly  affect. 

The  block  diagram  of  the  TKE  budget  (Figure  3.1)  is  presented  in  terms  of 
depth-averaged,  dimensional  variables  that  are  defined  in  terms  of  their  nondimensional 
counterparts  in  Table  6.  The  depth-averaged  variables  are  independent  of  mixed  layer 
depth,  unlike  the  nondimensional  variables.  Thus,  they  are  better  suited  for  comparing 
models  when  the  predicted  mixed  layer  depths  do  not  coincide.  The  wind  stress 
production  also  becomes  an  explicit  term.  These  characteristics  will  be  particularly 
useful  in  Chaper  IV  where  the  block  diagram  will  be  used  to  discuss  model  simulations 
at  OWS  P  and  OWS  N. 

The  TKE  budget  for  the  OPBL  is  shown  in  Figure  3.1.  The  terms  in  the  first  two 
equations  in  Table  3,  4  or  5  are  represented  by  the  arrows  in  the  figure.  The  third 
equation  in  all  three  tables  defines  P*.  Therefore,  no  new  processes  are  introduced  into 
the  TKE  budget  by  the  third  equation.  Wind  mixing  (G)  is  a  source  of  horizontal 
TKE.  Surface  cooling  (H  =  u*b*<0)  is  the  only  source  of  vertical  TKE.  The  sinks  of 


TABLE  6 

DIMENSIONAL,  DEPTH-AVERAGED  VARIABLES 


Total  TKE: 

E  = 

m3u» 

E*  h 

( 3. 10a) 

Vertical  TKE: 

W2  = 

h 

( 3. 10b) 

Surface  Buoyancy  Flux: 

H  = 

m,u*3 

H*— 2-— 

h 

( 3. 10c) 

Entrainment  Rate: 

P  = 

p  m3U*3 

*  h 

( 3. lOd) 

Rotation  Stress: 

R  = 

m3u*3 

h 

(3. lOe) 

Pressure  Redistribution: 

n  = 

p2(E  -  3W2) 

(  3.  lOf ) 

I/J  i  C  ^  y*  a  c  o  ^  >•  ^  1 1  ^  4*  i  ^  w  • 

c?  = 

m3u*3 

/  •?  1  Orr  ^ 

nlilU  QwLCSO  r  1  UUU^  L.XUU* 

h 

Dissipation: 

D  = 

E1/2(E  +  Z) 

( 3. lOh) 

with  Z  *  0  only  for  the  ZSTAR  model. 


TKE  are  dissipation,  surface  heating  and  entrainment  cooling.  The  last  two  convert 
the  TKE  into  potential  energy,  while  dissipation  converts  the  energy  into  heat.  The 
TKE  that  is  dissipated  is  effectively  lost  from  the  dynamical  system  of  the  ocean. 
However,  the  TKE  used  for  entrainment  continues  to  be  available  to  the  large-scale 
dynamics,  because  it  influences  the  potential  energy  of  the  system. 

For  wind-driven  regimes,  the  major  TKE  source  generates  horizontal  TKE.  Since 
dissipation  is  isotropic,  the  horizontal  TKE  is  dissipated  at  twice  the  rate  of  the  vertical 
TKE.  Entrainment  cooling  and  surface  heating  are  sinks  of  vertical  TKE  only.  Also 
entrainment  occurs  only  to  the  extent  that  there  is  vertical  TKE  to  transport  TKE  to 
the  entrainment  zone.  Since  the  major  sources  and  sinks  occur  in  different  TKE 
components,  the  transfer  of  energy  between  the  horizontal  and  vertical  components  is 
very  important  to  the  energetics  of  the  OPBL.  The  mixed  layer  dynamics  are 
significantly  affected  since  entrainment  is  controlled  by  the  vertical  TKE.  The 
pressure-rate  of  strain  terms  provide  a  mechanism  for  this  transfer  that  depends,  to  first 


order,  on  the  anisotropy  of  the  turbulence.  An  additional  mechanism  is  provided  by 
the  rotation  stress  interaction,  which  depends  on  the  orientation  of  the  planetary 
vorticity  in  relation  to  the  relative  vorticity  of  the  large-scale  turbulent  fluctuations. 

The  pressure-rate  of  strain  terms  (II)  will  transfer  energy  from  the  horizontal  to 
the  vertical  component  for  wind-dominated  regimes  and  from  vertical  to  horizontal  for 
buoyancy  driven  cases  (surface  cooling).  The  direction  of  transfer  for  the  rotation 
stress  terms  (R)  depends  on  the  wind  direction  and  not  on  the  anisotropy  of  the 
turbulence.  As  shown  in  Figure  3.1,  the  rotation  stress  transfers  energy  from  the 
horizontal  to  the  vertical  TKE  for  easterly  winds  and  from  the  vertical  to  the 
horizontal  for  westerlies.  On  diurnal  to  synoptic  time  scales,  either  wind  direction  is 
possible.  However,  on  seasonal  to  annual  time  scales  the  wind  direction  will  be  more 
persistent.  Thus,  the  possibility  exists  for  variations  in  the  ocean  mixed  layer  response 
to  atmospheric  cyclones  and  anticyclones  depending  on  location  relative  to  the 
direction  and  position  of  each  storm. 

In  the  ZSTAR  model,  the  dissipation  (D)  is  increased  by  the  addition  of  Z  (see 
3.10h),  which  is  the  dimensional  counterpart  to  Z*.  Thus,  Z  afFects  the  TKE  budget 
isotropically  and  is  always  a  sink  of  TKE.  It  reduces  the  entrainment  rate  by  reducing 
the  magnitude  of  which  is  the  total  TKE.  .  he  pressure-rate  of  strain  (II)  is  the  only 
TKE  conversion  mechanism  in  the  ZSTAR  model  and  it  always  reduces  the  anisotropy 
of  the  turbulence.  In  the  RSTAR  model,  D  is  smaller  since  Z  =  0.  Therefore,  the 
TKE  is  dissipated  more  slowly,  which  implies  that  more  TKE  is  available  for 
entrainment  than  in  the  ZSTAR  model.  The  rotation  stress  (R)  is  a  TKE  conversion 
mechanism  similar  to  II,  but  it  is  independent  of  the  TKE  distribution.  If  the  winds 
arc  easterlies  (westward  surface  stress),  R  converts  horizontal  TKE  to  vertical  TKE, 
which  further  increases  the  entrainment  rate.  For  westerlies,  eastward  surface  stress,  R 
converts  vertical  TKE  to  horizontal  TKE,  which  reduces  the  entrainment  rate.  For 
regions  with  westerly  winds,  the  hypothesis  is  that  R  may  reduce  or  eliminate  the  need 
for  Z  in  the  dissipation  parameterization  to  prevent  excessive  mixed  layer  deepening. 
As  is  clearly  demonstrated  in  Figure  3.1.  the  means  of  reducing  the  entrainment  rate  is 
very  different. 

D.  COMPARISON  OF  RSTAR  AND  ZSTAR  MODELS 

The  nondimensional  entrainment  rate  is  a  function  of  the  nondimensional  surface 
buoyancy  flux  (H*)  and  the  nondimensional  depth  (R*  in  the  RSTAR  model  or  Z*  in 
the  ZSTAR  model).  For  H*>0  (net  surface  heating),  the  buoyancy  flux  is  a  sink  of 
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TKE  (Figure  3.1  and,  for  example,  equation  3. Sc)  and  the  mixed  layer  is  in  a  forced 
convection  or  stable  regime.  In  this  regime,  an  equilibrium  depth  is  possible  for  which 
the  surface  buoyancy  damping  of  TKE  balances  the  surface  wind  stress  production. 
For  H*<0  (the  free  convection  or  unstable  regime),  the  surface  buoyancy  flux  is  a 
source  of  TKE.  Therefore,  an  equilibrium  depth  is  not  possible  without  some  sink  of 
TKE.  such  as  dissipation  or  rotation  stress  for  westerly  wind  regimes.  The  lack  of  such 
a  steady-state  solution  was  problematic  with  previous  mixed  layer  models  (Chapter  II). 
The  ocean  mixed  layer  depth  has  a  cyclical  steady-state  on  an  annual  time  scale. 
Unless  an  equilibrium  depth  is  possible  for  H*<0,  this  cyclical  steady-state  solution 
requires  an  annual  net  surface  heating  or  a  horizontal  divergence  in  the  upper  ocean  to 
generate  upwelling  at  the  base  of  the  mixed  layer. 

The  P*  =  0  curve  is  the  solution  for  mixed  layer  retreat.  Therefore,  it  gives  the 
equilibrium  depth  as  a  function  of  the  surface  fluxes.  For  the  RSTAR  or  the  ZSTAR 
models,  the  P*  =  0  curve  crosses  the  R*  or  Z*  axis.  Therefore,  a  cyclical  steady-state  is 
possible  in  both  models  without  a  net  surface  heat  flux  or  upwelling.  For  the  ZSTAR 
model,  no  equilibrium  depth  is  possible  for  H*<  O(-0.3).  For  RSTAR,  an  equilibrium 
depth  is  theoretically  possible  for  all  H*.  Thus,  the  RSTAR  model  has  the  potential  to 
maintain  a  cyclical  steady-  state  for  greater  wintertime  surface  cooling  rates  than  does 


the  ZSTAR  model. 


The  entrainment  rate  decreases  with  increasing  depth  for  both  models  (i.c. 
dP*,  dZ*  <0  and  5P*,'dR*<0).  For  the  same  surface  fluxes,  deep  mixed  layers  will 
deepen  less  rapidly  than  will  shallow  mixed  layers.  However,  the  decreasing  curvature 
of  the  P*  solution  (Figure  3.2a)  for  increasing  Z*  (i.e.  d2P+  dZ+2  <Q)  implies  that  the 
asymptotic  value  P«(Z*-»co)  is  greater  than  zero.  Thus,  mixing  is  still  theoretically 
possible  for  very  deep  mixed  layers.  Changes  in  Z*  also  are  more  effective  for  shallow 
mixed  layers  than  for  deep  mixed  layers  since  the  curvature  is  greater  for  shallow  mixed 
layers.  Conversely,  the  entrainment  rate  is  reduced  more  rapidly  with  increasing  depth 
by  R*  (Figure  3.2b)  since  the  curvature  of  the  P*  solution  increases  with  increasing  R* 
(i.e.  d2P*/dR*2>0).  This  implies  that  changes  in  R*  are  more  effective  for  deep 
mixed  layers  than  for  shallow  mixed  layers.  Also,  as  discussed  previously,  P*  must  be 
zero  for  some  value  of  R*>0  for  all  H*.  Thus,  the  RSTAR  model  can  not  deepen 


indefinitely. 


The  positive  curvature  of  the  P*  solution  with  respect  to  R*  also  means  that  P* 
approaches  an  asymptotic  value  as  R*-*-^.  This  means  that  the  entrainment  rate  is 


bounded  for  easterly  winds,  for  the  RSTAR  model.  For  easterly  wind*  ar.d  zero  or 
positive  surface  buoyancy  flux,  the  RSTAR  model  predicts  that  the  equilibrium  mixed 
layer  depth  given  by  in  (1.7)  is  negative.  Thus  in  the  RSTAR  model,  an 
equilibrium  mixed  layer  depth  can  not  be  achieved  in  regions  of  easterly  wind  without  a 
positive  surface  heat  flux  or  a  three-dimensional  balance  of  fluxes.  For  the  ZSTAR 
model,  only  the  condition  Z*£0  prevents  P,>-»30  for  negative  Z*. 

The  ratio  W*2,  E*  is  a  measure  of  the  anisotropy  of  the  integral  scale  TKE.  As 
discussed  previously,  the  large-scale  turbulence  can  not  be  isotropic  since  all  the 
production  and  loss  processes  are  not  isotropic  (see  Figure  3.1).  Only  the  dissipation 
affects  the  TKE  budget  isotropically.  If  the  large-scale  turbulence  were  isotropic. 
W'*2  E*  would  be  13.  For  wind-driven  regimes,  the  source  of  TKE  is  the  shearing  in 
the  mean  flow  generated  by  the  surface  wind  stress.  In  this  case,  TKE  is  produced 
only  in  the  horizontal  components  and  W*2  E*  should  be  less  than  1  3.  1  or 

buoyancy-driven  regimes,  the  TKE  production  is  also  anisotropic  since  the  source  of 
TKE  is  limited  to  the  vertical  TKE  component.  Also  the  entrainment  buoyancy  flux  is 
a  sink  of  vertical  TKE.  In  this  case,  W*2  E*  should  be  greater  than  1  3. 

The  ratio  W*2,E*  is  greater  for  the  ZSTAR  model  (Figure  3.3a)  than  for  the 
RSTAR  model  (Figure  3.3b)  for  all  Z*  and  for  R*  greater  than  zero.  This  is  consistent 
with  the  larger  values  of  P*  in  the  ZSTAR  model.  For  both  models,  W*2  E*  increases 
with  decreasing  H+  as  the  production  of  TKE  by  surface  buoyancy  fluxes  becomes  the 
primary  source  of  TKE  in  the  mixed  layer.  The  decrease  of  W*2  E*  with  increasing  Z* 
(Figure  3.3a)  is  much  slower  than  with  increasing  R*  (Figure  3.3b)  and  for  H*~-0.S 
W*2/E*  is  virtually  independent  of  Z*.  The  greater  variation  of  W*2  E*  with  R*  is 
because  R*  directly  alters  the  isotropy  of  the  TKE  distribution  (see  Table  3  or  Figure 
3.1)  by  transferring  TKE  between  the  eastward  and  vertical  components.  On  the  other 
hand,  Z*  reduces  the  TKE  by  enhancing  the  dissipation,  which  is  assumed  to  be 
isotropic.  Thus,  the  effect  of  Z*  on  P*  is  more  apparent  than  its  effects  on  W*2  E* 
since  reducing  E*  directly  reduces  P*  (3.8c).  However,  if  the  distribution  of  TKE  is  not 
initially  isotropic,  the  isotropic  dissipation  will  make  the  distribution  more  anisotropic 
with  time.  Thus,  for  W*2/E*<  1/3  (H* >  -0.8)  Z*  >0  tends  to  further  reduce  W*2  E*. 
The  trend  increases  with  increasing  Z*.  Thus  rate  of  decrease  of  W*2  E*  increases 
with  increasing  Z*  and  the  rate  of  decrease  of  W*2  E*  increases  as  W*2  E*  approaches 


zero. 


These  analyses  of  the  nondimensional  equations  for  the  TKE  models  show  that 
the  functional  form  of  the  entrainment  rate  (P*)  as  a  function  of  the  rotational 
dissipation  enhancement  (Z*)  is  very  different  than  as  a  function  of  rotation  stress 
(R*).  The  entrainment  rate  decreases  to  zero  with  increasing  rotation  stress  in  the 
RSTAR  model.  Rotation  stress  is  most  effective  in  altering  the  entrainment  rate  for 
deep  mixed  layers  having  large  rotation  stress  since  d2P*,5R*2>  0.  Conversely,  the 
entrainment  rate  in  the  ZSTAR  model  approaches  a  constant  with  increasing 
rotational  dissipation  enhancement.  Also,  the  dissipation  enhancement  more 
effectively  alters  P*  for  shallow  mixed  layers  with  small  rotational  dissipation 
enhancement  since  c2 P*,'<3Z*2  <  0. 

The  ratio  of  the  vertical  to  the  total  TKE  (W*2/E)  is  a  measure  of  the  isotropy  of 
the  TKE  distribution  among  the  components.  The  isotropy  of  the  TKE  is  affected 
directly  by  rotation  stress  through  the  transfer  of  TKE  between  the  east-west  and 
vertical  TKE  components.  However,  the  dissipation  alters  the  TKE  distribution  only 
to  the  extent  that  the  degree  of  isotropy  of  the  dissipation  rate  and  the  TKE 
distribution  do  not  match.  If  the  dissipation  removes  TKE  isotropically  from  the 
components  and  if  the  TKE  is  not  isotropically  distributed,  the  TKE  distribution  will 
become  more  anisotropic  with  time.  Thus,  the  ratio  W+2/E  was  found  to  vary  more  as 
a  function  of  rotation  stress  than  as  a  function  of  rotational  dissipation  enhancement. 


IV.  SIMULATIONS 


The  models  developed  in  Chapter  III  are  solved  numerically  for  E*.  \V*2  and  P* 
given  the  surface  fluxes  of  momentum  and  buoyancy,  and  for  the  RSTAR  model  the 
wind  direction  at  hourly  time  intervals.  During  entrainment,  the  model  density  profile 
is  mixed  to  a  depth  h  at  which  the  potential  energy  gained  by  entraining  denser, 
nonturbulent  water  from  below  the  OPBL  equals  the  vertical  convergence  of  TKE  from 
the  OPBL  to  the  entrainment  zone.  This  depth  need  not  correspond  to  a  level  nAz  of 
the  model,  where  n  is  an  integer.  Thus,  these  are  layer  rather  than  level  models.  Tor 
retreat,  the  mixed  layer  depth  is  set  to  an  equilibrium  depth  based  on  the  surface 
fiuxes.  Then  the  profile  between  that  depth  and  the  previous  depth  is  adjusted  at  each 
intervening  level  to  conserve  heat,  buoyancy  and  potential  energy.  Thus,  the  density 
profile  in  the  model  consists  of  a  well-mixed  layer  of  depth  h  below  which  the  density 
can  vary  as  a  function  of  z  at  intervals  of  Az  =  1  m.  The  maximum  depth  for  the 
profile  is  200  m. 

The  long  time  series  of  OWS  observations,  up  to  thirty  years  at  some  locations, 
allow  multiple  annual  simulations.  Shorter  case  studies  also  can  be  repeated  in  various 
seasons  and  years  for  intercomparison.  The  disadvantage  is  that  the  data  are  relatively 
coarsely  sampled  and  crudely  measured  when  compared  with  the  best  measurements 
that  can  be  taken  today.  Because  these  higher  quality  data  sets  result  from  intense 
field  experiments  that  last  only  weeks  to  months,  important  events  can  be  missed  and 
there  is  no  information  on  seasonal  variability.  For  these  reasons  and  due  to  its 
accessibility,  the  OWS  data  were  chosen  for  this  research  because  it  provides  the 
greatest  number  and  largest  variety  of  events  for  study. 

Both  the  weather  stations  used  for  these  studies  are  at  midlatitudes,  although 
equatorial  locations  would  have  been  preferable.  The  only  time  series  at  an  equatorial 
location  that  was  available  for  this  study  is  from  Gan  Island  (1°S,  73°E)  in  the  Indian 
ocean.  The  data  set  is  only  29  months  long  and  contains  several  data  gaps.  The  island 
is  located  in  the  climatological  path  of  the  westward  flowing  North  Equatorial  Current 
in  the  winter  and  the  eastward  flowing  Summer  Monsoon  Current  in  the  summer 
(Tchernia,  1980).  There  may  be  significant  three-dimensional  effects  with  large 
seasonal  variability  in  the  data.  In  fact,  in  the  Gan  Island  time  series  there  is  a  strong 


2  cycle  per  year  variation  in  the  zonal  wind  field  which  drives  a  zonal  jet  in  the  upper 
ocean  as  indicated  by  the  2  cycle  per  year  variation  in  the  upper  ocean  currents  which 
is  in  phase  with  the  wind  variations  (McPhadden,  1982). 

Bulk  formulas  are  used  to  calculate  the  surface  fluxes  from  the  meteorological 
observations.  The  values  of  the  coefficients  in  these  formulas  are  known  to  van'  with 
latitude  and  the  variations  are  not  as  well  documented  at  low  latitudes  as  at 
midlatitudes  {Reed,  1985).  The  meteorological  observations  were  taken  from  a  weather 
station  on  the  island  and  the  BT  casts  were  taken  several  miles  alternately  north  or 
south  of  the  island  at  weekly  or  less  frequent  intervals  (McPhadden,  1982).  Thus,  the 
meteorology  may  suffer  from  topographical  effects  due  to  the  island  and  the  BT  casts 
were  not  collocated  with  the  meteorological  observations.  For  these  reasons,  the  Gan 
Island  data  set  was  rejected  for  this  study. 

OWS  P  (50  °\,  145°  W)  and  OWS  N  (30°N,  140°  W)  provide  time  series  of  data 
needed  for  mixed  layer  simulations  on  annual  time  scales.  The  time  series  selected  for 
OWS  P  is  fifteen  years  long  and  for  OWS  X  eighteen  months.  For  the  simulations  in 
this  study,  the  surface  fluxes  of  heat  and  momentum  are  calculated  from  three  hourly 
OWS  observations  and  interpolated  to  hourly  intervals.  Since  evaporation  and 
precipitation  are  not  measured  at  ocean  weather  stations,  the  surface  salinity  fiux  will 
be  neglected.  Thus,  the  surface  buoyancy  flux  is  calculated  from  the  net  surface  heat 
flux.  The  density  as  a  function  of  depth  for  the  upper  200  m  of  the  ocean  is  needed  to 
initialize  and  verify  the  simulations.  The  OWS  observations  include  BT  casts,  which 
provide  ocean  thermal  structure  at  5  m  intervals  throughout  the  upper  ocean  at  twice 
daily,  or  more  frequent,  intervals.  These  can  be  used  to  calculate  density  profiles  for 
model  initialization  and  verification  if  the  salinity  is  known.  The  initial  profiles  for  the 
annual  simulations  are  composed  of  a  temperature  profile  representative  of  the  average 
thermal  conditions  for  early  January1  for  each  year  and  a  climatological  salinity  profile. 
For  OWS  P,  the  salinity  profile  contains  a  halocline  between  the  mixed  layer  and 
150  m.  The  profile  above  the  mixed  layer  depth  (based  on  the  temperature  profile)  has 
a  constant  salinity  of  34°/#0,  and  the  salinity  below  150  m  is  constant  at  35°, '00.  For 
OWS  N,  the  entire  profile  is  isohaline  with  a  salinity  of  35°  00.  Because  verification 
profiles  are  calculated  from  BT  casts  only  the  temperature  portion  of  the  density  profile 
is  verified. 

Recall  that  ZSTAR  is  the  Garwood  (1977)  model,  which  includes  the  rotational 
dissipation  enhancement  (Z*).  RSTAR  is  the  rotation  stress  (R*)  model  with  no 


dissipation  enhancement,  i.e.  Z*  =  0.  HSTAR  is  the  model  with  no  planetar;  rotation 
processes,  that  is  R#  =  Z*  =  0.  A  set  of  annual  simulations  will  be  made  using  the  same 
empirical  tuning  constants  in  all  three  models  (RSTAR,  ZSTAR  AND  HSTAR).  In 
this  tvay,  the  relative  responses  of  the  models  to  realistic  combinations  of  surface 
forcing  and  ocean  thermal  structure  can  be  compared.  Differences  in  the  responses  of 
the  models  can  be  attributed  directly  to  the  differences  in  the  modeled  physics. 

Several  periods  of  approximately  seven  days  in  duration  are  selected  for  further 
study  from  these  annual  simulations.  The  hypothesis  that  R*  can  cause  significant 
variations  in  diurnal  to  synoptic  scale  shallowing  and  deepening  events  will  be  tested  by 
comparing  the  mixed  layer  depths  predicted  by  the  RSTAR,  ZSTAR  and  HSTAR 
models  for  these  periods.  The  hypothesis  that  the  isotropy  of  the  TKE  distribution  for 
the  RSTAR  model  is  different  from  that  of  the  ZSTAR  model  will  also  be  tested  using 
these  synoptic  periods.  For  this  test,  the  terms  in  the  TKE  budget,  as  detailed  in 
Figure  3.1,  will  be  computed.  Another  working  hypothesis  is  that  the  rotation  stress, 
R*,  in  the  RSTAR  model  may  reduce  or  eliminate  the  need  for  Z*,  the  rotational 
dissipation  in  the  ZSTAR  model.  This  hypothesis  will  be  tested  by  tuning  seperatelv 
the  model  constants  to  optimize  simulations  of  the  annual  cycle  of  mixed  layer  depth 
and  temperature  for  the  RSTAR  model  and  the  ZSTAR  model.  Finally,  monthly 
statistics  and  annual  plots  of  mixed  layer  depth  and  temperature  will  be  used  to 
compare  the  seasonal  to  annual  time  scale  differences  between  the  RSTAR  and 
ZSTAR  models. 


A.  SURFACE  FLUX  CALCULATIONS 


The  surface  forcing  for  these  studies  was  calculated  using  the  formulas  and 
constants  discussed  in  Martin  (1985).  Some  of  the  surface  forcing  was  obtained  in  flux 
form  from  Mr.  Paul  Martin.  The  use  of  the  same  bulk  formulas  for  the  rest  of  the 
forcing  provided  a  consistent  data  set  for  the  simulations.  This  will  also  allow  further 
comparison  with  Martin's  results  in  the  future. 

The  wrind  stress  was  calculated  using  the  drag  coefficient  of  Garratt  (1977)  with 
the  wind  speed  corrected  to  a  10-m  height  using  a  logarithmic  profile  of  turbulent 
fluxes  in  the  atmospheric  surface  layer  after  the  method  of  Large  and  Pond  (1981). 
The  clear  sky  insolation  was  calculated  using  the  Milankovitch  formula  (List,  1984) 
with  the  Reed  (1977)  cloud  cover  correction  for  OWS  N  and  the  Tabata  (1964) 
correction  for  OWS  P.  The  use  of  two  different  cloud  cover  formulas  accounted  for  the 
variation  in  the  average  cloud  type  and  height  between  OWS  P  and  OWS  N  (Paul 


56 


W'wtJ 


j 


Martin,  personal  communication).  Net  longwave  radiation  was  calculated  using  the 
Beriiand  formula  (Wyrtki,  1965).  Latent  and  sensible  heat  fluxes  were  calculated  with 
bulk  formulas  using  the  exchange  coefficients  of  0.00122  at  OWS  X  and  0.00149  at 
OWS  P  as  was  done  by  Martin  (1985).  Solar  radiation  absorption  was  calculated  using 
the  parameterization  of  Garwood  (1976)  in  which  50%  of  the  irradiance  is  absorbed  in 
the  first  meter,  and  the  remainder  is  absorbed  exponentially  with  depth.  The  extinction 
depth  was  22  m  for  OWS  N  (corresponding  to  type  I  water,  Jerlov,  1976),  and  12.5  m 
for  OWS  P  (Jerlov's  type  II  water). 

The  bulk  formulas  and  the  values  of  the  exchange  coefficients  chosen  for  the 
surface  heat  and  momentum  flux  calculations  affect  the  performance  of  the  mixed  layer 
models.  The  choice  of  particular  formulas  and  coefficients  for  computing  the  boundary 
conditions  is,  in  essence,  part  of  the  tuning  process  for  the  models.  In  this  study,  the 
tuning  involved  constants  that  are  associated  with  the  parameterizations  in  the  ocean 
models,  rather  than  by  adjusting  the  atmospheric  forcing  formulas.  This  allows  the 
tuning  to  be  done  with  respect  to  the  physical  processes  that  are  parameterized  in  the 
models.  The  models  can  then  be  compared  using  the  same  surface  fluxes  and  the 
differences  in  model  simulations  can  be  ascribed  to  differences  in  model  physics. 

Martin  (1985)  performed  sensitivity  studies  of  the  bulk  formulas  and  exchange 
coefficients  for  OWS  P  and  OWS  N.  Gallacher  et  al.  (1983)  performed  sensitivity 
studies  for  the  solar  radiation  absorption  formulation.  Gallacher  ei  al.  (1983)  also 
compared  this  parameterization  with  several  others,  including  the  arctangent  model  of 
Zaneveld  and  Spinrad  (1980)  and  the  multiple  wavelength  model  of  Kondo  (1979),  and 
concluded  that  these  more  involved  parameterizations  did  not  significantly  alter  the 
depth-dependent  ocean  heating  due  to  shortwave  radiation.  The  choices  of  formulas 
and  coefficients  for  this  study  were  based  on  those  results.  The  choices  provide  a 
model-independent  set  of  surface  fluxes  with  the  minimum  annual  bias  in  net  surface 
heating  that  is  consistent  with  direct  observations  of  surface  fluxes. 

B.  DESCRIPTION  OF  OCEAN  WEATHER  STATION  DATA 

Annual  simulations  were  conducted  for  OWS  P  (50°N,  145°W)  and  OWS  N 
(30°N,  140°W).  At  OWS  P,  the  years  1961,  1965,  1966  and  1967  were  chosen  because 
these  years  each  had  more  than  2000  BT  casts  that  could  be  used  for  model 
verification.  In  1965,  1966  and  1967,  the  observations  wrere  distributed  rather 
uniformly  throughout  the  year  (Figures  4. lb-4. Id).  In  1961  (Figure  4.1a).  the  BT 
observations  were  clustered  around  several  intensive  observation  periods  for  the 
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Internal  Wave  Experiment  (IWEX).  Furthermore,  1961  was  the  year  used  by  Martin 
(1985)  for  an  intercomparison  of  the  ZSTAR  model,  the  Niiler  (1975)  model  and  the 
Mellor  and  Yamada  (1974)  level  2  and  2  1/2  models,  at  OWS  P  and  OWS  \.  The 
ZSTAR  model  predicted  the  seasonal  cycle  better  than  did  the  Mellor  models  or  the 
Niiler  model.  The  ZSTAR  model  provided  a  good  simulation  at  OWS  N  but  tended  to 
shallow  too  much  in  late  summer  and/or  deepen  too  rapidly  in  the  fall  at  OWS  P, 
depending  on  the  model  constants.  Thus,  the  ZSTAR  model  provides  a  good 
comparison  model  since  it  is  as  good  or  better  than  other  state  of  the  art  vertically- 
integrated  models.  The  observations  for  1961  at  OWS  N  are  distributed  uniformly 
throughout  the  year  (Figure  4.1e).  However,  there  are  only  1650  observations 
compared  to  an  average  of  2350  for  the  years  at  OWS  P. 

OWS  P  is  in  a  region  of  predominately  westerly  winds.  OWS  X  is  in  a  region  of 
predominately  easterly  winds  during  the  spring,  summer  and  fall.  The  wind  direction  at 
OWS  X  is  more  variable  during  the  winter  w’hen  the  subtropical  high  pressure  center 
moves  south  to  approximately  30°X,  140°W  from  its  summertime  postion  of 
approximately  40°X,  150°W  (Tchemia,  1980).  For  1961,  the  wind  direction  at  OWS 
X  is  southerly  to  southwesterly  in  January  and  early  February  then  predominately 
easterly  for  the  rest  of  the  year.  These  stations  are  separated  by  20°  latitude.  Thus,  it 
might  seem  they  provide  an  appropriate  test  for  the  wind  direction  and  latitude 
dependence  of  rotation  stress.  However,  the  conditions  at  OWS  P  and  OWS  X  differ 
in  more  than  these  aspects  (Elsberry  and  Camp,  1978).  The  magnitude  of  the  wind 
forcing  and  surface  heating  varies  considerably  from  OWS  P  to  OWS  X  as  shown  in 
Table  7.  The  values  in  Table  7  suggest  that  the  heat  budget  for  1961  at  OWS  X  may 
be  more  one-dimensional  than  at  OWS  P,  since  the  magnitude  of  the  annual  net 
surface  heating  (17  W  m'2)  at  OWS  X  is  less  than  the  magnitude  of  the  annual  net 
surface  cooling  (25  W  m‘2)  at  OWS  P  for  the  same  year.  There  is  no  evidence  at  either 
location  of  a  climatological  heating  or  cooling  trend.  Therefore,  any  net  heating  or 
cooling  from  the  atmosphere  must  be  balanced  by  a  horizontal  or  vertical  heat  flux  in 
the  ocean. 

The  climatological  oceanographic  conditions  indicate  that  OWS  P  is  located  in  a 
region  more  likely  to  be  dominated  by  the  local  atmospheric  forcing  than  is  OWS  X. 
Because  OWS  X  is  located  within  the  region  of  influence  of  the  subtropical  front, 
horizontal  advection,  convergence  and  divergence  are  expected  to  be  more  significant 
(Tchernia,  1980).  Indeed,  there  are  indications  in  the  temperature  time  series  at  OWS 
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N  that  the  subtropical  front  oscillates  across  the  station.  Climatologically,  the  ocean 
surface  currents  at  OWS  X  are  from  the  north  and  northeast.  This  current,  which  is 
part  of  the  California  Current  system,  transports  colder  water  into  the  region  of  OWS 
X.  Thus,  a  net  cooling  due  to  horizontal  transport  is  possible.  At  OWS  P,  the 
climatological  currents,  which  form  the  northern  edge  of  the  subtropical  gyre,  are  from 
the  west.  The  currents  are  weaker  on  average  than  at  OWS  X  and  the  transport  is 
along  the  isotherms  rather  than  nearly  norma!  to  the  isotherms  as  at  OWS  X. 
Therefore,  the  net  heat  transport  is  expected  to  be  smaller  at  OWS  P  than  at  OWS  X 
on  average.  Thus,  the  wind  direction  and  latitude  are  only  two  of  many  differences 
between  these  ocean  weather  stations.  Differences  between  the  simulations  and 
observations  will  reflect  all  the  variations  in  atmospheric  forcing  and  oceanic 
conditions.  However,  by  choosing  periods  of  diurnal  to  synoptic  time  scales,  during  late 
winter  to  early  spring  the  effects  of  planetary  rotation  are  likely  to  be  most  apparent. 

Martin  (1985)  observed  that  the  ZSTAR  model  predicted  a  mixed  layer  depth  at 
OWS  P  that  was  slightly  deeper  than  observed  during  the  summer  and  early  fall, 
whereas  the  mixed  layer  depth  was  well  predicted  at  OWS  X.  The  equilibrium  mixed 
layer  depth  for  the  RSTAR  model  is  given  by 
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and  for  the  ZSTAR  model  by 
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The  ratio  of  rotation  stress  to  surface  buoyancy  flux  is  given  by  (1.9) 

<i>  -  avT»  -  Jl; 

"  pE7(0)  "  H*  ' 

In  Garwood  et  al.  (1985a),  we  showed  that  typical  summertime  values  are  <I>  =  0.19  at 
OWS  P  and  <l>  =  0.08  at  OWS  X.  These  values  have  the  correct  magnitude  and  sign  to 
improve  the  mixed  layer  depth  prediction  at  both  OWS  P  and  OWS  X  because  the 
resulting  equilibrium  depth  for  the  RSTAR  model  relative  to  that  for  the  ZSTAR 
model  is  smaller  at  OWS  P  than  at  OWS  X.  The  value  of  R*  (3.6f)  relative  to  the 
surface  buoyancy  flux  H*  C?.6e),  rather  than  the  absolute  magnitude  of  R*.  determines 
the  efficacy  of  R*  for  altering  the  equilibrium  depth. 
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TABLE  7 

ANNUAL  AVER.AGE  SURFACE  FLUXES  AT  OWS  P  AND  X 


Tx  Ty 

Qs 

Qb 

Qe 

Qh 

Qn 

OWS  P 
1961 

(X  in2) 
0.102  0.013 

216. 

(W  m' 

94. 

2) 

86. 

n. 

-25. 

1965 

0.099  0.013 

221. 

98. 

103. 

19. 

-1. 

1966 

0.0S0  -0.005 

216. 

95. 

99. 

19. 

-3. 

1967 

0.100  0.012 

223. 

99. 

108. 

17. 

1. 

OS  N 
1961 

-0.030  -0.006 

349. 

111. 

227. 

28. 

17. 

t  is  the  east-west  component  of  surface  wind  stress 
t  is  the  north-south  component  of  surface  wind  stress 
Qs  is  the  solar  radiation  (positive  downward) 

Qb  is  the  longwave  radiation  (positive  upward) 

Qe  is  the  evaporative  heat  flux  (positive  upward) 

Qh  is  the  sensible  heat  flux  (positive  upward) 

Qn  is  the  net  heat  flux  (Qb  +  Qe  +  Qh-Qs) 

C.  SYNOPTIC-SCALE  INTERCOMPARISONS  OF  MODELS 

Several  synoptic-scale  periods  are  chosen  from  annual  simulations  of  1965  and 
1967  at  OWS  P  and  from  the  annual  simulations  of  1961  at  OWS  N  to  demonstrate 
the  latitudinal  and  seasonal  variations  among  the  RSTAR,  ZSTAR  and  HSTAR 
models  and  the  time  scales  of  these  variations.  As  discussed  in  the  introduction,  R* 
decreases  with  increasing  latitude.  Conversely,  Z*  increases  with  increasing  latitude. 
The  time  scale  of  variations  in  the  rotation  stress  ranges  from  the  inertial  scale  of  the 
planetary  rotation  to  synoptic  scale  of  the  variations  in  wind  direction. 

Both  Z*  and  R*  are  proportional  to  h.  Because  the  mixed  layer  is  at  or  near  a 
deep  equilibrium  state  during  late  winter,  diurnal  to  synoptic  scale  events  that  occur 
during  late  winter  or  during  the  spring  transition  are  most  likely  to  be  affected  by  R* 
and  Z*.  The  spring  transition  is  that  period  when  the  mixed  layer  changes  from  the 
deep  equilibrium  state  to  the  shallow  summertime  regime.  This  change  can  occur  quite 
abruptly  due  to  the  modulation  of  the  mixed  layer  depth  by  the  diurnal  cycle  of  surface 
heating.  The  transition  occurs  on  a  day  when  the  wind  stress  production  of  TKE  is 
not  sufficient  to  return  the  mixed  layer  to  the  deep  equilibrium  depth  from  which  the 
diurnal  retreat  began.  The  date  of  this  transition  is  significantly  correlated  at  OWS  P 
w’ith  the  accumulation  of  heat  above  the  seasonal  thermocline  (Elsberry  and  Garwood. 


1978).  Rotation  stress  may  affect  this  transition  by  augmenting  (easterly  wind  regimes) 
or  reducing  (westerly  wind  regimes)  the  entrainment  rate.  The  frequency  of  diurnal  and 
synoptic-scale  cycles  of  shallowing  and  deepening  in  the  RSTAR  simulations  are  quite 
different  from  those  for  the  ZSTAR  model  during  the  late  winter  to  early  spring  for  all 
five  years  of  this  study  (Figures  4.10a  through  4.14a). 

In  the  summer,  the  winds  are  light,  the  surface  heating  is  large  and  the  insolation 
is  near  the  annual  maximum.  Consequently,  the  layer  is  probably  too  shallow  for  R* 
or  Z*  to  significantly  affect  the  TKE  budget.  During  the  fall  deepening  period,  wind 
mixing  and  surface  cooling  dominate  the  TKE  budget  and  the  entrainment  rate  is 
primarily  determined  by  the  net  TKE  production  and  the  density  gradient  below  the 
mixed  layer.  The  direct  effects  of  R*  are  expected  to  be  dominated  by  the  TKE 
production  due  to  wnnd  stress  and  surface  cooling,  and  the  TKE  loss  due  to 
entrainment  cooling  and  dissipation.  Although  variations  in  the  density  profile  of  the 
seasonal  pycnocline  caused  by  rotation  stress  during  the  spring  transition  could  later 
affect  the  fall  deepening,  demonstrating  such  seasonal  teleconnections  would  be  very 
difficult.  For  these  reasons,  the  winter- spring  period  was  examined  for  examples  of  the 
influences  of  rotation  stress.  Synoptic  periods  characteristic  of  the  summer  and  fall 
conditions,  as  discussed  above,  were  also  chosen  to  demonstrate  the  seasonal 
variability. 

For  these  simulations,  the  values  of  the  tuning  constants  are  a  compromise 
among  the  values  that  were  determined  for  the  RSTAR  and  ZSTAR  models  by 
optimizing  the  simulations  of  the  annual  cycle  at  OWS  P  for  1965,  and  the  values 
determined  from  theoretical  estimates  and  from  independent  experimental  results  (see 
the  next  section  for  details  concerning  the  tuning  and  the  optimal  values  for  each 
model).  The  value  for  Pj,  which  is  a  measure  of  entrainment  efficiency,  is  0.3.  It 
controls  the  size  of  the  entrainment  cooling  relative  to  the  dissipation  (Figure  3.1). 
The  strength  of  the  pressure  redistribution  term  relative  to  the  dissipation  is  p2,  which 
is  equal  to  0.5  for  these  simulations.  The  relative  efficiency  of  TKE  production  due  to 
wind  and  breaking  surface  waves  is  m3  (m3  =  7.5).  The  rotation  stress  term  does  not 
contain  a  tuning  constant,  as  can  be  seen  in  Table  2.  The  coefficient  for  Z*  is  p3.  It  is 
measure  of  the  dissipation  enhancement  due  to  rotation  (Table  2).  For  these 
simulations,  p3  is  set  equal  to  1.0.  The  depth-averaged,  dimensional  variables  defined 
in  equations  (3.10a)-(3.10f)  of  Table  2  will  be  used  for  the  comparisons  in  this  section. 
These  variables  differ  from  the  nondimensional  variables  in  that  they  are  independent 
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of  mixed  layer  depth  and  the  wind  stress  production  becomes  an  explicit  term. 

(G  =  mjU*3). 

An  example  of  the  winter  through  spring  transition  period  is  the  seven  days  taken 
from  the  annual  simulations  of  1967  at  OWS  P  (Figure  4.2).  A  most  striking  feature  of 
this  period  is  that  the  RSTAR  and  HSTAR  models  deepen  much  more  rapidly  than  the 
ZSTAR  model  (Figure  4.2a).  After  100  days  of  simulations,  all  three  models  predict 
similar  mixed  layer  depths.  Then  the  ZSTAR  model  deepens  to  90  m  on  day  103.5 
w'hile  the  HSTAR  and  RSTAR  models  deepen  through  the  permanent  pycnocline,  to  a 
depth  of  175  m  in  the  same  time  span.  This  divergence  of  the  entrainment  rates  (P)  is 
not  due  to  variations  in  the  simulated  density  profiles  below'  the  mixed  layer,  as  will  be 
demonstrated  below.  It  results  from  the  reduced  dissipation  (D)  in  the  HSTAR  and 
RSTAR  models  relative  to  that  in  the  ZSTAR  model  (Figure  4.2b).  The  TKE  that  is 
made  available  by  this  reduction  in  D  is  used  to  increase  P  (Figure  4.2c).  The  result  is 
that  the  mixed  layer  in  the  ZSTAR  model  deepens  more  rapidly,  while  the  TKE 
remains  approximately  the  same  for  all  three  models  (Figure  4.2g). 

An  alternative  hypothesis  is  that  the  structure  of  the  seasonal  pycnocline 
immediately  below  the  mixed  layer  in  the  HSTAR  and  RSTAR  simulations  is  different 
from  that  of  the  ZSTAR  simulation  due  to  variations  in  mixed  layer  depth  and  density 
during  the  preceding  100  days  of  the  simulations.  If  that  is  the  case,  the  rates  of  mixed 
layer  deepening  and  the  ultimate  equilibrium  depths  are  different  due  to  the  differences 
in  the  pycnocline  structure  into  which  the  mixed  layer  is  entraining  rather  than  to 
variations  in  the  mixed  layer  TKE  budget.  To  test  this  hypothesis,  the  RSTAR  and 
HSTAR  simulations  were  recalculated  starting  on  day  101  with  the  profile  for  that  day 
and  time  generated  by  the  ZSTAR  model.  The  results  (not  shown)  are  similar  to 
Figure  4.2.  The  HSTAR  and  RSTAR  models  deepened  more  than  the  ZSTAR  model 
and  achieved  a  final  equilibrium  depth  of  approximately  150  m.  This  is  shallower  than 
the  175  m  equilibrium  depth  in  Figure  4.2a.  However,  the  variations  in  deepening  rate, 
dissipation  rate,  vertical  and  total  TKE  are  sufficiently  similar  to  those  in  Figure  4.2  to 
conclude  that  they  are  primarily  the  result  of  the  variations  in  the  TKE  budget,  rather 
than  variations  in  the  underlying  thermocline. 

The  crosses  in  Figure  4.2a  (and  in  the  subsequent  mixed  layer  depth  plots)  are 
the  mixed  layer  depths  calculated  from  the  OWS  data.  For  the  data,  the  mixed  layer 
depth  is  defined  as  the  greatest  depth  with  a  temperature  less  than  0.I°C  colder  than 
the  sea-surface  temperature.  The  trend  of  the  mixed  layer  depth  predicted  by  the 


ZSTAR  model  agrees  with  the  data  better  than  does  that  predicted  by  the  HSTAR  and 
RSTAR  models  (Figure  4.2a).  However,  the  most  interesting  times  for  comparisons 
among  the  RSTAR,  ZSTAR  and  HSTAR  models  are  the  diurnal  retreat-deepening 
events  (Figure  4.2a)  around  days  104,  105,  106  and  107.  During  these  events,  the 
mixed  layer  is  quite  deep  compared  to  typical  values  for  the  summer  or  fall  (Figures  4.5 
and  4.6)  even  at  the  shallowest  point  of  these  retreat  events  (approximately  100  m  for 
the  RSTAR  and  HSTAR  models).  Thus,  the  rotation  stress  and  the  rotationally 
enhanced  dissipation  will  have  maximum  effect  on  the  dynamics  of  the  mixed  layer. 
The  diurnal  events  predicted  by  all  three  models  are  not  clearly  represented  in  the  data. 
There  does  appear  to  be  some  diurnal  period  variations  in  the  data,  particularly  around 
days  105  and  106,  which  have  a  larger  amplitude  than  that  predicted  by  the  ZSTAR 
model.  The  amplitude  of  the  variations,  roughly  50  to  100m,  appears  to  agree  better 
with  the  amplitude  of  the  events  predicted  by  the  HSTAR  and  RSTAR  models  than 
with  that  predicted  by  the  ZSTAR  model.  The  ZSTAR  model  appears  to  fit  the  long¬ 
term  trend  in  the  data.  However,  the  RSTAR  model  may  better  predict  the  diurnal 
events  that  appear  in  the  data,  although  the  average  mixed  layer  depth  in  the  RSTAR 
model  is  clearly  too  deep. 

There  is  a  large  diurnal  to  synoptic  time  scale  signal  in  W2/E  (Figure  4.2h)  since 
W’  =  0  whenever  the  mixed  layer  retreats.  In  the  RSTAR  and  HSTAR  models,  this 
ratio  is  two  to  six  larger  than  in  the  ZSTAR  model.  However,  this  increase  is  not 
accompanied  by  an  increase  in  the  total  transfer  rate,  Fl-R,  (Figure  4.2i)  as  might  be 
expected  from  Figure  3.1.  In  fact,  the  total  transfer  rate  for  the  RSTAR  and  HSTAR 
models  is  less  than  for  the  ZSTAR  model.  Thus,  the  increase  in  W2E  is  not  the  result 
of  increased  transfer  of  TKE  from  the  horizontal  to  the  vertical  component.  It  is  a 
direct  result  of  the  decreased  dissipation  (see  below). 

Since  the  dissipation  is  specified  to  be  isotropic,  dissipation  of  the  vertical  TKE 
accounts  for  1,3  of  the  total  TKE  dissipated  (Figure  3.1).  The  vertical  TKE  is 
approximately  5%  of  the  total  TKE  in  the  HSTAR  and  RSTAR  models  and  only  1% 
in  the  ZSTAR  model  (Figure  4.2h).  Therefore,  removing  a  third  of  the  total  dissipated 
TKE  from  the  vertical  component  may  account  for  the  difference  between  the 
estimates  from  the  experimental  results  and  the  model  results.  If  this  is  the  case,  the 
TKE  distribution  may  be  less  isotropic  in  the  simulations,  particularly  in  the  ZSTAR 
simulation,  than  it  is  in  nature.  This  tentative  conclusion  can  not  be  verified  with  the 
data  being  used  in  this  study.  Thus  it  must  be  left  as  a  hypothesis  to  be  tested  in 
future  work. 
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The  diurnal  mixing  events  can  be  explained  in  terms  of  the  IKE  budget  shown  in 
Figure  3.1.  The  mixed  layer  depth  is  well  correlated  with  H,  the  dimensional,  depth- 
averaged  surface  buoyancy  flux  (3.10c).  The  shallowing  corresponds  to  the  increasingly 
positive  H  that  results  from  the  diurnal  increase  in  insolation.  The  maximum  H  occurs 
at  the  minimum  depth  of  the  cycle,  which  corresponds  to  W2,  E  =  0,  and  H  acts  as  a 
sink  of  vertical  TKE.  Although  the  existence  and  vertical  range  of  the  events  do 
depend  on  G  (3.10g),  the  depth  of  the  mixed  layer  is  poorly  correlated  with  G  and  D 
throughout  the  events.  The  average  value  of  G  is  approximately  2.0xl0'3cmV3.  This 
is  the  major  source  of  TKE  during  most  of  the  events.  With  the  primary  source  of 
TKE  in  the  horizontal  component  and  the  major  sink  in  the  vertical,  the  total  TKE 
transfer  (IT  -  R)  achieves  a  local  maximum  value  of  1.0-2.0xl0*3cm2s‘J  during  each 
event.  In  all  cases,  P  =  0  until  the  minimum  depth  is  acliieved  and  then  P  increases  to 
a  maximum  value  of  approximately  1.5xl0‘4cm2s'3  at  the  point  of  maximum  depth. 
The  behavior  of  P  mimics  that  of  W2/E,  as  expected.  The  maximum  value  of  W2  E  is 
about  0.1  to  0.2.  For  the  ratio  of  rotation  stress  to  surface  buoyancy  flux  (O>0)  the 
RSTAR  model  should  retreat  to  a  shallower  depth  than  the  HSTAR  model,  for  which 
R  =  0.  This  is  confirmed  in  Figure  4.2a,  since  the  mixed  layer  depth  for  the  RSTAR 
model  is  always  shallower  than  for  HSTAR. 

When  variations  in  G  and  R  are  in  phase,  the  wind  direction  must  be  easterly. 
This  is  generally  the  case  in  the  events  discussed  here.  At  the  start  of  the  event  on  day- 
104.  R  (Figure  4.2e)  and  G  (Figure  4.2d)  are  no  larger  than  in  other  events.  However, 
the  event  corresponds  to  a  local  maximum  for  R  and  local  minimum  for  G.  This 
implies  that  the  wind  speed  is  reduced  and  the  direction  shifted  such  that  the  winds 
were  blowing  from  the  east.  All  the  other  events  are  local  minimums  for  both  R  and 
G.  This  extra  conversion  of  vertical  to  horizontal  TKE  due  to  the  increasing  R  during 
day  104  is  coupled  with  a  reduction  of  TKE  due  to  decreased  G  throughout  the  period. 
Consequently,  the  RSTAR  model  retreats  30  m  whereas  the  mixed  layer  depth  in  the 
HSTAR  model  does  not  change.  Conversely,  both  models  retreat  equally  for  the  event 
of  day  107  because  R  is  approximately  zero  throughout  the  period  (Figure  4.2e).  The 
smaller  rate  of  deepening  of  the  RSTAR  model  compared  to  the  HSTAR  model  for 
days  101  through  103  is  a  direct  result  of  the  reduced  total  transfer,  IT  -  R,  (Figure 
4.2i)  due  to  the  positive  R  (Figure  4.2e).  This  difference  in  the  rate  of  deepening  also 
appeared  in  the  test  for  which  the  RSTAR  and  HSTAR  models  were  restarted  at  day- 
101.  Therefore,  it  is  not  a  consequence  of  any  differences  in  the  pycnocline  structure 
below  the  mixed  layer. 


Another  case  during  the  winter  through  spring  transition  is  from  OWS  P  during 
days  72-80  of  1965  (Figure  4.3).  The  most  significant  difference  between  this  case  and 
the  previous  one  (Figure  4.2)  is  that  R  is  negative  throughout  most  of  the  period 
because  the  winds  are  easterly.  As  a  result,  the  mixed  layer  depth  (Figure  4.3a)  from 
the  RSTAR  simulation  is  deeper  than  that  from  HSTAR  by  approximately  10%.  The 
ZSTAR  model  is  shallower  than  either  the  HSTAR  or  RSTAR  model  and  the  mixed 
layer  depths  calculated  from  the  BT  casts  are  generally  closer  to  the  ZSTAR  solution. 
However,  there  still  appears  to  be  more  diurnal  to  synoptic  scale  variability  in  the  data 
than  in  the  ZSTAR  model  predicted  depths.  For  ZSTAR,  the  average  W2  E  (Figure 
4.3h)  is  small,  about  0.01,  as  in  the  previous  case.  For  HSTAR,  W2/E  is  about  twice 
the  previous  value.  For  the  RSTAR  model,  W2,  E  is  20%  larger  than  for  the  HSTAR 
model.  This  results  from  the  enhanced  transfer  of  TKE  from  the  horizontal  to  the 
vertical  component.  In  this  case,  R*  acts  to  increase  the  isotropy  of  the  TKE  budget. 
During  the  synoptic  period  from  day  74  through  day  77,  the  mixed  layer  from  the 
RSTAR  simulation  becomes  increasingly  deeper  than  the  mixed  layer  depth  from 
HSTAR  (Figure  4.3a)  since  the  winds  are  easterly  throughout  the  period.  As  the 
mixed  layer  deepens  and  the  entrainment  rate  shows  the  average  W2  E  (Figure  4.3h)  is 
decreases. 

A  final  example  during  the  winter-spring  transition  is  from  OWS  N  during  days 
77-84  of  1961  (Figure  4.4).  For  this  case,  R  oscillates  between  positive  and  negative 
values  with  a  two  to  three  day  period  (Figure  4.4e).  The  positive  peaks  in  R 
correspond  to  peaks  in  G  (Figure  4.4d);  whereas  R<0  events  corresponds  to  low  wind 
speeds.  For  this  case,  the  mixed  layer  depths  (Figure  4.4a)  predicted  by  all  three 
models  are  similar,  due  to  the  relatively  large  production  of  TKE  by  the  surface 
buoyancy  flux.  In  this  example,  H  (Figure  4.40  is  more  negative  by  approximately  a 
factor  of  two  than  for  the  1967  OWS  P  case  (Figure  4.20-  The  mixed  layer  depths 
calculated  from  the  BT  casts  do  not  show  the  diurnal  cycles  predicted  by  all  three 
models.  No  model  fits  the  data  appreciably  better  than  the  other  two.  The  mixed 
layer  from  the  RSTAR  simulation  is  shallower  than  from  the  HSTAR  simulation  when 
R  is  positive  on  day  78  and  deeper  on  day  79  when  R  is  negative.  The  ratio  W2  E 
(Figure  4.4h)  is  greater  than  0.33  for  ail  the  deepening  events  which  indicates  that  the 
events  are  buoyancy-driven.  For  the  deepening  event  on  day  79.  W2  E  from  the 
RSTAR  model  is  larger  than  for  the  HSTAR  or  ZSTAR  model.  This  is  a  direct  result 
of  the  enhanced  transfer  of  TKE  from  the  horizontal  to  the  vertical  component  due  to 
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the  negative  R.  The  isotropy  of  both  the  RSTAR  and  HSTAR  models  arc  greater  than 
for  the  ZSTAR  model  due  to  the  reduced  dissipation  in  those  models.  At  the  start  of 
day  79,  the  dissipation  (Figure  4.4b)  in  the  ZSTAR  model  is  roughly  50%  larger  than 
in  the  RSTAR  and  HSTAR  models. 

It  was  hypothesized  that  the  effects  of  R  and  Z  would  be  the  greatest  during  the 
winter  through  the  spring  transition  period  discussed  aboved.  The  effects  of  planetary 
rotation  should  be  smaller  during  summer  since  the  mixed  layer  is  relatively  shallow 
and  the  winds  are  light.  A  synoptic  period  from  days  210  to  217  at  OWS  P  during 
1967  will  be  used  to  verify  this  hypothesis.  For  this  period  (Figure  4.5).  the  three 
models  give  virtually  identical  results.  The  average  for  R  during  this  period  is  almost 
zero  and  the  maximum  value  of  3.0xl0'5cm2s‘3  (Figure  4.5e)  is  more  than  an  factor  of 
10  smaller  than  the  maximum  value  for  the  winter  case  (Figure  4.2e).  The  mixed  layer 
depth,  0(3  m)  (Figure  4.5a).  and  the  average  ratio  of  vertical  to  total  TKE  of  0.07 
(Figure  4.5h)  are  virtually  identical  for  all  three  models. 

In  the  fall  (Figure  4.6),  R  and  Z  are  more  effective  than  in  the  summer. 
However,  the  mixed  layer  deepening  is  accelerating  mainly  in  response  to  increasingly 
large  wind  stress  and  positive  surface  buoyancy  flux.  For  the  period  in  Figure  4.6,  days 
250  to  257  from  OWS  P  for  1967.  the  average  wind  stress  production  is 
G  =  1.0xl0‘2cm2s'3  (Figure  4.6g)  and  the  surface  buoyancy  flux  (Figure  4.60  is 
-6.0xl0'4cm2s'J  £  H  SS  4.0xl0"4cm2s"3,  with  negative  values  indicating  production  of 
TKE.  By  contrast,  the  average  R  is  only  5.0xl0'5cm2s'3.  During  this  period,  all  three 
models  are  deepening  at  similar  rates.  The  rate  of  mixed  layer  deepening  is  limited  not 
only  by  the  rate  of  production  of  TKE  by  the  surface  fluxes  but  also  by  the  strength  of 
the  pycnocline.  For  example,  the  fairly  substantial  storm  on  day  255  (maximum 
G=  3.0xl0'2cm2s'3  and  a  maximum  H  *  -4.5xl0'4cm2s*3)  only  deepened  the  mixed 
layer  5m.  The  mixed  layer  depth  from  the  ZSTAR  model  is  approximately  25% 
shallower  than  that  of  the  HSTAR  or  RSTAR  models  throughout  the  period  and 
W2  E  is  smaller  because  of  the  enhanced  dissipation  in  the  ZSTAR  model.  However, 
the  rate  of  deepening  0(1  m,  day)  is  approximately  the  same  for  all  three  models  and 
the  difference  between  the  RSTAR  and  HSTAR  models  is  negligible  even  for  the  large 
negative  R  (-4.0xl0'4cm2s'3)  on  day  254.5.  This  indicates  that  the  rate  of  transfer  of 
TKE  between  the  horizontal  and  vertical  components  is  not  a  significant  part  of  the 
TKE  budget  during  this  period.  These  synoptic  periods  are  typical  of  the  seasons 
represented. 
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To  summarize  the  results  of  the  cases  discussed  above,  both  R*  and  Z*  were 
shown  to  be  relatively  ineffectual  during  the  summer  and  fall.  During  the  winter 
through  the  spring  transition.  R*  and  Z*  substantially  alter  the  diurnal  to  synoptic 
scale  mixing  events  and  the  long  term  trend  in  mixed  layer  depth.  Variations  in  the 
dynamics  of  the  upper  ocean  can  be  described  in  terms  of  the  TKE  budget  of  the 
OPBL  (Figure  3.1).  The  isotropy  of  the  TKE  distribution  for  the  OPBL  changes 
substantially  on  diurnal  to  synoptic  time  scales.  The  turbulence  is  likely  to  be  most 
anisotropic  for  retreat  events  (W*2->0).  The  isotropy  of  the  TKE  budget  predicted  by 
the  RSTAR  model  is  considerably  different  than  that  predicted  by  the  ZSTAR  model, 
as  is  the  entrainment  rate.  For  the  RSTAR  model,  the  isotropy  (e.g.,  Figure  4.2h)  is 
larger  than  for  the  ZSTAR  model  by  as  much  as  a  factor  of  10  and  the  entrainment 
rate  (e.g.,  Figure  4.2c)  is  2  to  2.5  times  greater  than  for  the  ZSTAR  model.  These 
results  from  the  synoptic  scale  studies  appear  to  contradict  the  nondimensional  results. 
The  iatter  predicted  greater  variation  in  the  isotropy  (Figure  3.3b)  and  reduced 
entrainment  rate  (Figure  3.2b)  as  a  function  of  R*,  particularly  for  deep  mixed  layers. 
This  apparent  contradiction  is  due  to  the  smallness  of  Z*  and,  particularly,  R*  in  the 
synoptic  studies.  For  the  winter-spring  transition  cases  most  values  were  in  the  range 
O^Z*  ^0.5  and  0£R*^0.1.  In  the  nondimensional  plots  the  values  are 
0^Z*,R*:S2.0  for  P*  and  0:SZ*^1.0,  -0.5  ^R*^  0.5  for  W*2/E.  The  maximum 
values  used  in  the  nondimensional  plots  are  well  within  the  range  of  reasonable  values 
at  OWS  P  and  OWS  N  but  these  larger  values  do  not  occur  for  the  near  equilibrium 
situations  shown  in  the  synoptic  scale  figures. 

At  least  for  the  parameter  values  used  in  the  synoptic  studies,  R*  did  not  provide 
sufficient  damping  of  P*  to  prevent  excessive  mixed  layer  deepening  for  the  winter 
cases.  However,  the  ZSTAR  model  predicts  a  significantly  less  isotropic  TKE  balance 
in  the  OPBL.  The  isotropy  of  shear-forced  mixing  layers  can  be  estimated  from 
various  laboratory  studies  that  have  been  reported  in  the  literature  (Table  8).  There 
are  no  field  measurements  in  the  upper  ocean,  from  which  estimates  of  the  isotropy  of 
the  integral  scale  turbulence  can  be  made.  Based  on  the  estimates  in  Table  S.  the 
average  value  of  W*2/E'v0.05  for  the  RSTAR  and  HSTAR  models  probably  is  an 
underestimate  and  the  value  of  W*2  E  ~0.01  for  the  ZSTAR  model  very  likely  is  an 
underestimate  ( Figure  4.2h).  However,  the  laboratory  measurements  are  for  wall- 
bounded,  pure  shear  flows.  The  OPBL  is  bounded  by  a  free  surface  and  there  can  be 
significant  production  and  damping  of  the  TKE  by  buoyancy  fluxes.  For  the  case  of 


buoyant  damping,  the  isotropy  of  the  OPBL  should  be  less  than  the  estimates  from  the 
laboratory  measurements,  which  are  for  neutrally  buoyant  flows. 

TABLE  8 

ESTIMATES  OF  THE  ISOTROPY  OF  WALL-BOUNDED  SHEAR  LAYERS 


Source 

Type  of  Flow 

W.2 

Smits  et  al.  (1979) 

Convex  curvature 

0.16 

Muck  et  al.  (1985) 

Convex  curvature 

0.22 

Townsend  (1976) 

Flat  plate 

0.25 

Hunt  and  Joubert  (1979) 

Convex  curvature 

0.25 

So  and  Mellor  (1973) 

Convex  curvature 

0.29 

The  temporal  and  spatial  resolution  and  the  precision  of  the  OWS  data  sets  used 
in  this  study  was  not  sufficient  to  verify  the  increased  frequency  and  magnitude  of 
diurnal  to  synoptic  scale  shaliowing-deepening  events  predicted  by  the  RSTAR  model. 
Further  simulations  arc  recommended  using  data  from  experiments  such  as  the  Storm 
Transfer  and  Response  Experiment  (STREX)  and  the  Mixed  Layer  Dynamics 
Experiment  (MILDEX)  to  verify  the  effects  of  rotation  stress.  These  experiments  were 
intense  synoptic  period  observations  of  the  upper  ocean  and  the  surface  fluxes  at 
higher  temporal  resolution  and  with  greater  accuracy-  and  precision  than  the  OWS 
data.  However,  these  data  sets  may  not  contain  events  with  the  east-west  wind  stress, 
which  is  needed  to  test  the  variations  among  the  models.  A  set  of  experiments  a  ned 
specifically  at  verifying  the  rotation  stress  effects  may  be  required  Measuremn  of 
surface  stress  to  ±  l.OxlO"2  X  m'2  and  surface  heat  fluxes  to  ±  10  W  m'*  are  now 
possible  (Large  tt  al.,  1986).  These  should  produce  sufficiently  accurate  surface  fluxes 
to  test  the  models  introduced  in  this  research.  The  capability  now  exists  to  measure 
temperature  profiles  to  millidefree  precision.  This  would  allow  the  differences  in 
simulated  mixed  layer  temperature  and  depth  among  the  models  to  be  tested  along 
with  the  assumption  of  a  "well-mixed  similarity  profile. 

The  magnitude  of  the  processes  that  affect  the  TKE  budget  (e  g  R.  D.  and  fl* 
calculated  for  the  synoptic  scale  studies  vary  from  IO*5cmV’  to  lO*2cmV\  and  the 
magnitude  of  the  TKE  is  0(  IO‘2cm2s'2)  The  difference  between  the  values  of  W,’  f  , 
predicted  by  the  RSTAR  and  ZSTAR  models  implies  a  difference  in  vertical  IKE 
averaged  over  the  mixed  layer  of  2xlO*JcmV2  Differences  in  the  magnitudes  ol  R.  D. 


II  and  the  other  processes  shown  in  Figure  3.1  range  from  10*5cm2s'3  to  10‘3cmV3 
between  the  RSTAR  and  ZSTAR  models.  The  lower  limit  for  measurements  of  the 
rate  of  dissipation  of  TICE  in  the  ocean  is  2.0xl0*6cm2s'3  (Osborn  and  Lueck,  1984) 
due  to  instrument  noise.  Thus,  the  differences  in  the  TKE  budgets  calculated  in  this 
work  may  be  measureable  with  current  technology.  These  dissipation  measurements 
resulted  from  microstructure  measurements  of  turbulent  velocity  shear,  absolute  values 
of  turbulent  velocity  may  not  be  measurable  to  the  same  precision.  Also  horizontal 
transects  through  the  mixed  layer  w'ould  be  required  to  collect  time  series  of  sufficient 
duration  to  obtain  stable  correlations  of  the  turbulent  velocity  components.  Such 
transects  are  difficult  even  using  a  relatively  stable  underwater  platform  such  as  a 
submarine  (Rolf  Lueck,  personal  communication). 

These  studies  of  synoptic  periods  selected  from  annual  simulations  with  the 
RSTAR,  ZSTAR  and  HSTAR  models  demonstrated  that  for  the  same  values  of  the 
tuning  constants  the  models  predicted  quite  different  responses  to  the  same  forcing. 
Rotation  stress  changes  the  frequency  and  magnitude  of  diurnal  to  synoptic  scale 
shallowing  and  deepening  events  predicted  for  deep  mixed  layers.  During  these  events, 
the  depth  changes  calculated  by  the  RSTAR  model  were  10%  to  30%  greater  than 
those  calculated  by  the  HSTAR  model.  Simulations  from  both  these  models  were 
substantially  different  from  simulations  with  the  ZSTAR  model  during  the  winter, 
during  spring  transition  and  dunng  the  fall  deepening.  All  three  model  calculated 
similar  summertime  mixed  layers.  The  TKE  distributions  predicted  by  the  RSTAR  and 
HSTAR  models  were  10%  to  40%  larger  and  more  isotropic  than  that  predicted  by  the 
ZSTAR  model.  This  was  the  result  of  the  larger  dissipation  associated  with  the 
inclusion  of  the  rotational  dissipation  enhancement  in  the  ZSTAR  model.  The 
dissipation  is  assumed  to  be  isotropic  in  all  three  models,  whereas  W*2  E*  was  usually 
less  than  0. 1  in  the  simulations.  Laboratory  measurement  indicate  that  the  values  for 
W*:  L*  are  probably  too  small. 

D.  TUNING  USING  ANNUAL  SIMULATIONS 

For  the  simulations  discussed  in  the  previous  section,  the  values  of  the  tuning 
constants  were  based  primarily  on  theoretical  estimates  and  experimental  results.  As 
discussed,  the  ZSTAR  model  simulations  best  lit  the  trends  in  the  data.  In  general, 
any  of  the  models  could  be  tuned  to  estimate  the  trend  in  the  data  lor  synoptic-scale 
periods  This  is  due  partis  to  the  relatives  narrow  range  ol  surface  (luxes  and  mixed 
liver  variations  that  occur  over  such  short  periods  However,  annual  periods  include  a 


sufficiently  broad  range  of  physical  processes  that  tuning  alone  is  not  sufficient  if  the 
model  physics  is  inadequate  or  incomplete.  As  discussed  in  Chapter  II,  the  need  for  a 
mechanism  to  prevent  excessive  wintertime  mixing  in  bulk  models  is  well  known.  The 
purpose  of  this  and  the  following  section  is  to  demonstrate  the  requirement  for  such  a 
mechanism  and  to  determine  the  extent  to  which  rotation  stress  may  provide  this 
mechanism  for  regions  of  westerly  wind.  The  results  of  Garwood  et  al.  (1985b)  suggest 
the  possibility  that  rotation  stress  may  eliminate  the  need  for  rotationally  enhanced 
dissipation  on  seasonal  to  annual  time  scales.  The  hypothesis  for  this  section  is  that, 
for  the  proper  values  of  the  tuning  parameters,  the  RSTAR  model  can  simulate  the 
annual  cycle  of  the  OPBL  for  westerly  wind  regimes  without  excessive  winter  deepening 
or  a  generally  poor  fit  to  the  annual  cycle. 

The  RSTAR  and  ZSTAR  mixed  layer  models  are  optimally  tuned  for  this  part  of 
the  study.  That  is,  the  model  constants  are  adjusted  in  an  attempt  to  obtain  a  best  fit 
to  the  observed  mixed  layer  depth  and  temperature.  A  more  complete  discussion  of  the 
definition  of  best  fit  between  the  models  and  observations  is  given  below.  The  rotation 
stress  term  does  not  contain  a  tuning  constant.  However,  altering  the  constants  in  the 
other  terms  implicitly  alters  the  effects  of  R*.  Because  the  models  are  nonlinear  the 
effects  of  the  physical  processes  are  not  independent  and  all  the  constants  must  be 
adjusted  if  the  models  are  to  be  optimally  tuned  for  a  given  set  of  physical  processes. 

The  parameters  that  will  be  tuned  and  the  processes  they  affect  are  shown  in 
Figure  3.1.  For  the  models  discussed  in  this  work,  the  efficiency  of  TKE  production 
due  to  wind  and  breaking  surface  waves  and  the  convergence  of  TKE  within  the  mixed 
layer  is  m3.  The  entrainment  efficiency  is  pj.  The  strength  of  the  pressure 
redistribution  term  is  determined  by  p2.  In  the  ZSTAR  model,  the  p3  coefficient  for  Z* 
is  a  measure  of  the  strength  of  the  dissipation  enhancement.  The  coefficients  p,,  p2 
and  p3  are  all  defined  relative  to  the  lowest  order  parameterization  of  dissipation, 
E*3  2.  The  RSTAR  model  contains  three  tuning  parameters  whereas  the  ZSTAR 
model  contains  four. 

The  1965  data  at  OWS  P  were  chosen  for  these  tuning  experiments  because  of 

the  large  number  of  evenly  distributed  observations  (Figure  4.1c).  Furthermore,  the 

one-dimensional  mixed  layer  model  assumptions  are  met  since  the  net  annual  heating 

was  relatively  small  during  that  year  (Table  7).  The  models  were  tuned  by  choosing  a 

* 

set  of  statistical  measures  to  be  minimized  on  an  annual  time  scale.  The  statistical 
parameters  chosen  were  the  bias, 


Xbias  =  (-Xm  *  “Xd)/N 
and  the  root  mean  square  error  (RMSE). 

Xrmse  =  {2C(  Xm  “  Xd>2/*Q}l/2 

of  the  mixed  layer  temperature  (T)  and  depth  (h).  In  the  above  definitions,  the  xm's 
are  the  simulated  values  and  the  xd's  are  the  observations.  The  summations  are  over 
all  BT  casts  for  the  year  and  the  concurrent  model  value  is  specified  at  the  hour  nearest 
the  time  of  the  BT  cast.  The  bias  measures  the  fit  of  the  phase  and  amplitude  of  the 
low  frequency  (seasonal  to  annual  time  scale)  variability  between  the  observations  and 
the  simulations.  Matching  the  timing  and  peak  values  of  the  summertime  maximum 
temperature  and  minimum  depth  and  the  wintertime  minimum  temperature  and 
maximum  depth  between  the  observations  and  the  simulations  will  minimize  the  bias. 
The  RMSE  is  minimized  by  matching  the  phase  and  amplitude  of  the  high  frequency 
(synoptic  to  diurnal  period)  variability  in  the  model  and  the  observations. 

For  a  specified  value  of  p2,  the  bias  in  the  RSTAR  model  can  be  minimized  for 
the  locus  of  (m3,  p,)  points  defining  the  zero  contour  lines  on  Figures  4.7a  and  4.8a  for 
the  depth  and  temperature,  respectively.  If  p2  is  increased,  the  zero  contour  is  shifted 
toward  smaller  (m3,  p}),  but  the  curvature  remains  roughly  the  same.  Thus,  the  zero 
bias  for  temperature  or  depth  forms  a  surface  in  (m3,  Pj,  p2)  space  that  decreases  with 
increasing  p2  and  has  concave  curvature  in  the  (m3,  pL)  plane.  If  p2  is  increased,  the 
transfer  of  TKE  from  the  horizontal  to  the  vertical  increases.  Thus,  the  same  mixed 
layer  depth  and  temperature  could  be  achieved  with  less  TKE  production  (m3)  and  less 
effective  mixing  (pj).  For  fixed  p2,  increasing  the  TKE  production  (m3)  implies  that 
less  efficient  mixing  (reduced  pt)  is  needed  to  achieve  the  same  mixed  layer  depth 
(MLD)  and  mixed  layer  temperature  (MLT). 

Although  the  temperature  (Figure  4.8b)  and  depth  (Figure  4.7b)  RMSE  s  have  a 
pattern  similar  to  that  of  the  biases,  they  do  not  overlay  the  zero  bias  contours. 
Rather  the  RMSE  minima  correspond  to  biases  for  which  the  model  is  warmer  (Figure 
4.8a)  and  shallower  (Figure  4.7a)  than  indicated  by  the  observations.  A  shallower 
OPBL  will  tend  to  have  more  high  frequency  variability  due  to  the  reduced  thermal 
inertia.  The  rather  large  RMSE's  indicate  that  the  model  variations  lack  some  of  the 
high  frequency  variability  that  exists  in  the  BT  data.  The  high  frequency  variations  in 
the  BT  data  may  be  due  to  such  factors  as  noise,  motion  of  the  ship  on  station, 
horizontal  advection  of  the  temperature  patchiness  that  has  been  observed  in  the  upper 
ocean  (R.  W.  Garwood,  personal  communication),  or  it  may  be  due  to  the  differences 
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in  definition  of  MLD  and  MLT  for  the  model  and  the  data.  Some  bias  in  model  versus 
data  comparisons  can  be  expected  from  the  difference  between  the  observational 
definition  of  MLD  and  MLT  given  above  and  the  model  MLD  and  MLT.  In 
particular,  the  large  minimum  RMSE  and  the  displacement  of  the  RMSE  minimum 
relative  to  the  zero  bias  is  consistent  with  this  difference  in  definition.  This  study  is 
concerned  more  with  comparing  the  physical  processes  parameterized  in  the  RSTAR 
and  ZSTAR  models  than  with  achieving  the  optimal  model  versus  data  comparison. 
Therefore,  no  attempt  will  be  made  to  correct  for  any  possible  biases  that  may  be 
associated  with  the  different  definitions  of  MLD  and  MLT. 

As  discussed  above,  the  minima  of  the  temperature  and  depth  biases  and 
RMSE's  for  the  RSTAR  model  form  surfaces  in  (m3  pp  p2)  space.  There  is  no 
indication  that  the  minimum  values  become  smaller  anywhere  in  parameter  space  or 
that  the  surfaces  converge.  Thus,  the  fit  between  the  RSTAR  model  and  the  data,  at 
least  based  on  these  statistics,  will  not  improve  significantly  for  some  set  of  parameter 
values  different  from  the  values  of  the  tuning  constants  used  in  the  previous  section. 
The  conclusion  is  either  that  R*  can  not  eliminate  Z*  at  midlatitudes  on  annual  time 
scales  or  that  these  annual  statistics  are  insufficient  to  allow  quantitative  tuning  of  the 
model. 

The  statistical  estimates  that  were  chosen  are  not  perfect  since  the  values  can  be 
minimized  for  the  wrong  reasons.  For  example,  the  values  of  the  biases  can  be 
minimized  by  a  model  annual  cycle  that  is  too  warm  in  the  summer  and  too  cold  in  the 
winter  such  that  the  average  of  the  biases  is  small.  The  RMSE  can  be  minimized  if  the 
dominant  time  scale  of  the  model  is  the  same  as  the  dominant  time  scale  of  the  data 
and  they  are  in  phase  even  if  the  amplitudes  are  significantly  different.  Thus,  it  is 
possible  to  achieve  a  “good”  fit  with  respect  to  these  statistics  for  a  model  that  is 
significantly  different  from  the  data.  These  statistics  do  not  appear  sufficient  to  form 
an  explicit,  quantitative  method  of  model  tuning.  Since  it  is  apparent  that  R«  is  not 
sufficient  to  eliminate  Z*,  no  attempt  will  be  made  to  determine  an  improved  set  of 
model  tuning  statistics  for  this  study.  However,  such  a  study  will  be  necessary  for 
future  work  involving  real  time  upper  ocean  predictability. 

Estimates  for  some  parameters  can  be  made  from  independent  experiments  and 
from  theory.  The  dimensional  equilibrium,  or  retreat,  depth  can  be  computed  for  the 
RSTAR  model  by  setting  P*-W«2-0  in  Table  3  and  substituting  for  the 
nor.dimensional  variables  from  Table  2  to  obtain  a  dimensional  retreat  depth.  The 
resulting  retreat  depth  is 


(p2  +  2/3  +  2<J>) 

which  is  the  same  as  (2.8)  in  Garwood  ct  al.  (1985a).  This  indicates  that  p2^  1.3  is 
required  for  hf^0.  This  establishes  a  lower  bound  on  p2,  at  least  for  the  RSTAR 
model.  Ze man  and  Tennekes  (1975)  summarized  laboratory  and  atmospheric 
observations  that  indicate  that  p2  should  be  approximately  12.  Garwood  ei  al. 
(1985a)  estimated  m3  to  be  approximately  seven  based  on  measuren  nts  by  Oakey  and 
Elliott  (1982).  Values  for  p3  and  Pj  are  less  quantifiable.  However,  p;  must  be  less 
than  one  and  p3  should  be  of  order  one.  If  p3  is  significantly  different  from  one.  it 
most  likely  represents  the  effects  of  physical  processes  other  than  dissipation,  and  these 
processes  should  be  modeled  explicitly.  The  model  does  not  seem  to  be  very  sensitive 
to  variations  of  pv  Preliminary  tuning  simulations  indicated  that  varying  p,  by  ±0.1 
altered  the  annual  average  depth  bias  approximately  ±2.5  m  about  an  annual  average 
bias  of  approximately  0  m  and  varied  the  annual  average  temperature  bias  roughly 
±0.0S°C  about  an  annual  average  average  bias  of  1.3 °C.  A  more  complete  study  of 
the  sensitivity  of  the  ZSTAR  model  to  variations  of  p3  might  indicate  what  physical 
processes  are  being  parameterized  with  Z*  and  how  the  parameterization  could  be 
improved. 

The  models  were  tuned  by  changing  m3  in  the  neighborhood  of  seven  and 
allowing  Pj  to  range  between  0.1  and  0.5  until  no  extreme  values  of  the  monthly  depth 
biases  and  RMSE  s  were  obtained  for  any  months  cf  1965  at  OWS  P.  The  resulting 
viIlcs  lor  the  model  constants  are  summarized  in  Table  9. 

TABLE  9 

MODEL  CONSTANTS  FOR  OWS  SIMULATIONS 


ZSTAR 


RSTAR 


E.  SEASONAL  AND  ANNUAL  SCALE  INTERCOMPARISONS  OF  MODELS 


The  problem  of  excessive  wintertime  deepening  in  early  bulk  models  of  the  OPBL 
was  discussed  in  Chapter  II.  Several  solutions  to  this  problem  were  reviewed.  In 
particular,  the  rotationallv  enhanced  dissipation  parameterization  of  Garwood  (1977) 
has  been  discussed  in  some  detail  and  is  a  major  component  of  the  ZSTAR  model 
being  used  in  this  study.  To  demonstrate  the  effectiveness  of  this  dissipation 
parameterization,  simulations  of  1965  at  OWS  P  with  the  HSTAR  and  ZSTAR  models 
are  compared  using  the  final  tuning  constant  values  chosen  for  the  ZSTAR  model 
(Table  9).  The  ZSTAR  and  HSTAR  simulations  are  similar  in  the  summer  (Figure 
4.9).  However,  the  HSTAR  simulation  is  significantly  deeper  than  the  ZSTAR 
simulation  or  the  observed  mixed  layer  for  days  0-130.  Also  the  fall  (days  275  to  325) 
deepening  rate  in  the  HSTAR  simulation  is  too  large.  This  results  in  a  wintertime 
mixed  layer  depth  that  is  approximately  75  m  deeper  than  the  observations  after  day- 
325.  Only  the  existence  of  the  halocline  between  the  initial  mixed  layer  depth  and 
150  m  prevents  the  mixed  layer  from  being  even  deeper  (Figure  4.9).  These  simulations 
clearly  demonstrate  the  need  for  some  additional  process  or  processes  to  reduce  the 
efficiency  of  mixed  layer  deepening  in  the  HSTAR  model  and  allow  a  cyclic  steady- 
state  for  the  OPBL  depth. 

The  seasonal  to  annual  time-scales  effects  of  rotation  stress  and  rotationallv 
enhanced  dissipation  are  in\estigated  using  annual  simulations  for  1961.  1965.  19<>6 
and  1967  a;  OWS  P  and  1961  at  OWS  N  The  RSTAR  and  ZSTAR  models  use  the 
values  of  the  tuning  constants  given  :n  Table  9.  The  biases  and  RMSE  s  are  calculated 
every  36  days  and  annually  for  all  the  simulations.  The  annual  cycics  of  predicted 
3-hourly  MLD  (Figures  4.10a  -  4.14a)  and  MLT  (Figures  4.10b  -  4.14b)  for  each  model 
arc  similar  in  gross  features  and  comparable  to  the  observed  values.  The  monthly 
biases  and  RMSE's  for  the  MLD  and  MLT  (Table  10  -  14)  also  have  some  common 
characteristics  on  annual  time  scales.  Both  models  tend  to  be  too  warm  and  shallow  in 
the  summer  compared  to  the  observations. 

Spring  transition,  which  is  the  usually  abrupt  transition  from  a  deep  wintertime 
regime  to  a  shallow  summertime  regime,  occurred  around  day  9<i  for  OWS  P.  1961. 
I9o6  and  1967.  about  day  70  tor  OW  S  P.  1965  and  OWS  N.  1901.  and  near  day  60  lor 
OWS  P.  1966.  The  models  tend  to  be  deeper  and  colder  than  the  observations  prior  to 
the  start  of  the  spring  transition.  This  may  be  due  in  part  to  the  lack  of  surface 
sal.nity  flux  in  the  boundary  conditions.  The  annual  evaporation  and  precipitation 
may  play  a  significant  role  in  limiting  the  wintertime  deepening  at  OWS  P  both  late  in 


the  year  and  prior  to  spring  transition  (Garzon,  1987).  Also,  the  lack  of  accurate 
temperature  and  salinity  profiles  for  initialization  is  important  to  model  performance 
prior  to  spring  transition.  Before  the  construction  of  this  seasonal  thermoclinc,  the 
model  MLD  is  extremely  sensitive  to  the  strength  of  the  pycnocline  immediately  below 
the  mixed  layer.  Thus,  the  initial  pycnocline  is  critical  to  the  model  performance 
during  this  period.  After  the  spring  transition,  the  model  begins  to  build  a  seasonal 
thermocline. 

Tables  10  through  14  give  the  monthly  and  the  annually  averaged  temperature 
and  depth  biases  for  the  ZSTAR  and  RSTAR  models.  For  OWS  P  during  1967,  the 
annual  average  surface  heat  flux  of  +0.8  W  m*2  (cooling)  is  the  smallest  of  the  five 
years.  For  both  models  (Table  10),  the  annual  temperature  biases  and  RMSE's  are 
smaller  than  for  any  other  year  at  OWS  P.  Both  models  have  annual  average 
temperature  biases  of  less  than  1.0°C.  The  ZSTAR  bias  was  0.4°C  (positive  bias 
means  the  model  was  warmer  than  the  observations),  and  the  RMSE  was  0.6°C.  For 
RSTAR,  the  bias  was  0.9°C  and  the  RMSE  was  1.0°C.  The  values  for  RSTAR  are 
about  twice  the  ZSTAR  values.  However,  the  annual-average  depth  bias  is 
significantly  smaller  for  RSTAR  and  the  annual-average  depth  RMSE's  are 
comparable.  The  depth  bias  was  -2.3  m  and  the  RMSE  was  33  m  for  RSTAR;  for 
ZSTAR,  the  bias  was  -15  m.  the  RMSE  was  40  m.  A  negative  bias  means  the  model 
was  shallower  than  the  observations.  The  ZSTAR  model  is  warmer  and  shallower  than 
the  observations  every  month  except  October  and  November,  which  accounts  for  the 
larger  depth  bias.  The  RSTAR  model  achieves  a  better  annual  average  bias  b\  having 
a  monthly  average  temperature  that  is  too  cold  in  the  winter  and  fall  and  too  warm  in 
the  summer.  The  maximum  monthly  temperature  bias  in  the  RSTAR  model  of  3.0 °C 
in  August  is  three  times  that  of  the  ZSTAR  model.  Similarly,  the  RSTAR  model 
mixed  layer  is  too  deep  in  the  winter  (maximum  bias  =32  m)  and  too  shallow  in  the 
summer  (maximum  bias  =  -6  m). 

The  RMSE's  for  the  RSTAR  and  ZSTAR  models  are  comparable  for  all  months 
except  January,  February  and  December  1967  (Table  10).  In  those  three  months,  the 
ZSTAR  model  has  significantly  more  diurnal  and  synoptic  scale  variability  than  either 
the  data  or  the  RSTAR  model  (Figure  4.13a).  The  diurnal  to  synoptic  variability  is 
similar  for  both  models  for  days  120  to  240.  For  days  60  to  120  and  240  to  330,  the 
RSTAR  mode!  has  less  diurnal  to  synoptic  variability  than  the  observations,  whereas 
the  ZSTAR  model  has  more.  Thus,  the  "comparable"  RMSE  values  for  the  two 
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models  are  achieved  quite  differently.  For  RSTAR,  the  variability  is  less  and  for 
ZSTAR  it  is  greater  than  in  the  observations. 

The  net  annual  surface  heat  flux  is  approximately  zero  in  all  three  years, 
1965-1967  at  OWS  P.  The  annual  temperature  bias  and  RMSE  for  the  RSTAR  model 
are  similar  for  all  three  years.  For  the  ZSTAR  model,  the  annua!  temperature  bias  and 
RMSE  is  larger  for  1965  and  1966  than  for  1967.  In  1965,  the  increase  is  due  mainly 
to  larger  values  early  in  the  year,  whereas  in  1966  the  increase  is  the  result  of  larger 
differences  between  the  model  and  the  observations  late  in  the  year.  Thus,  it  is 
possible  to  obtain  similar  values  for  the  annual  biases  and  RMSE's  with  quite  different 
monthly  results. 

For  1961  at  both  ocean  weather  stations,  the  annual  average  surface  heat  flux  is 
larger,  by  an  order  of  magnitude,  than  in  the  other  three  years.  At  OWS  P,  there  was 
an  annual  average  surface  warming  of  25  W  m"2.  At  OWS  N,  there  is  an  annual 
average  surface  cooling  of  17  W  m-2.  The  annual  average  temperature  and  depth 
biases  for  both  the  RSTAR  and  ZSTAR  simulations  are  two  to  three  times  larger  at 
OWS  P  for  1961  than  they  were  for  the  other  three  years.  Both  models  are  shallower 
and  warmer  than  the  observations.  For  ZSTAR,  the  annual  average  temperature  bias 
is  1.8°C,  which  is  three  times  the  values  for  1965  through  1967.  The  depth  bias  is 
-9.4  m,  which  is  similar  to  the  other  years.  Both  the  annual  average  temperature  bias 
for  RSTAR  of  3.2°C  and  the  depth  bias  of  -7.5  m  are  larger  than  the  values  from 
other  years.  These  results  are  consistent  with  the  large  annual  average  surface 
warming  for  1961  at  OWS  P. 

The  biases  and  RMSE's  for  OWS  N  are  almost  identical  for  both  models  and  are 
smaller  than  at  OWS  P  for  1961  and  comparable  to  OWS  P  for  1965-1967.  The  annual 
average  temperature  bias  is  -0.5  °C  for  both  models  and  the  model  depth  is  5  m 
shallower  than  the  observations  for  ZSTAR  and  4.3  m  shallower  for  RSTAR.  The 
annual  average  surface  stress  at  OWS  N  was  approximately  0.03  N  m"2,  compared  with 
an  average  of  about  0.1  N  m'2  for  the  four  years  at  OWS  P  (Table  7).  At  OWS  N,  the 
annual  average  north-south  stress  is  a  factor  of  five  smaller  than  the  east-west  stress. 
Thus,  there  is  a  relatively  large  annual  average  surface  cooling  with  light  easterly  winds 
at  OWS  N  compared  to  an  approximately  zero  net  surface  heating  and  moderate  to 
strong  westerlies  at  OWS  P.  This  implies  that  R*  and  Z*  will  be  smaller  at  OWS  N  due 
to  the  small  surface  stress  and  due  to  the  relatively  shallow  mixed  layer.  The  winter 
maximum  mixed  layer  is  0(120  m)  at  OWS  N  compared  to  (150  m)  at  OWS  P.  Also 


H*  will  be  relatively  large  and  negative  at  OWS  Thus,  the  relative  effects  of  Z*  and 
R*  at  OWS  N  will  be  small  and  mixed  layer  deepening  will  mainly  be  driven  by  the 
surface  buoyancy  flux.  The  similar  results  for  the  two  models  indicates  that  the  annual 
simulations  are  not  very  sensitive  to  Z*  and  R*  at  OWS  \  for  1961. 

The  result  that  the  simulated  depths  are  too  shallow  on  average  may  indicate 
that  mixing  due  to  surface  buoyancy  flux  is  underestimated  in  the  models.  The  same 
result  could  be  due  to  tuning  the  dissipation  to  be  too  large  to  compensate  for  too 
much  mixing  by  the  surface  stress.  The  simulations  are  colder  than  the  data  due  to  the 
large  surface  cooling  being  distributed  over  a  too  shallow  layer.  The  effects  due  to 
latitudinal  variations  in  R*  are  obscurred  by  the  variations  caused  by  the  significantly 
different  surface  buoyancy  and  momentum  fluxes  at  the  two  ocean  weather  stations. 

The  RSTAR  model  is  generally  warmer  and  shallower  in  the  summer  and  colder 
and  deeper  in  the  winter  compared  to  the  ZSTAR  model  and  the  observations  at  OWS 
P.  This  is  characteristic  of  an  incomplete  parameterization  of  dissipation  (Gaspar, 
1986).  At  the  same  time,  the  ZSTAR  model  has  more  diurnal  to  synoptic  scale 
variability  than  the  RSTAR  model  or  the  observations  (for  these  values  of  the  tuning 
constants.  Table  9).  This  can  be  the  result  of  too  much  dissipation  such  that  the  mixed 
layer  retreats  toe  much  for  slight  variations  in  wind  speed.  In  general,  the  ZSTAR 
model  simulated  the  data  better  than  the  RSTAR  model.  However,  the  ZSTAR  model 
rarely  achieved  a  good  fit  to  the  data  for  all  seasons  of  the  year.  Thus,  a  better  annual 
fit  had  to  be  achieved  at  the  expense  of  one  or  more  seasonal  periods. 

Based  on  the  results  of  the  synoptic  studies,  the  results  from  tuning  the 
parameters,  and  the  results  discussed  above,  the  conclusion  is  that  a  reduction  of 
mixing  efficiency  is  needed,  at  least  for  midlatitude  mixed  layers.  The  rotation  stress  is 
not  sufficient  to  prevent  excessive  mixed  layer  deepening  in  regions  of  westerly  wind  at 
midlatitudes.  However,  the  results  of  Garwood  ei  al.  (1985b)  indicate  that  R*  may  be 
sufficient  to  prevent  excessive  deepening  in  equatorial  regions  where  R*  approaches  a 
maximum  and  Z*  approaches  zero.  However,  over  much  of  the  tropics  the  winds  arc 
easterly  and  R*  would  enhance  mixing.  The  synoptic  studies  discussed  above  indicate 
that  for  westerly  wind  regions  R*  can  reduce  the  wintertime  mixed  layer  depth  10%  to 
15%.  Thus,  a  smaller  value  of  Z*  may  be  required  in  a  model  which  includes  both  the 
rotation  stress  and  the  rotationally  enhanced  dissipation.  A  smaller  value  of  Z*  may 
reduce  some  of  the  problems  in  the  ZSTAR  simulations  that  were  discussed  above.  An 
alternative  solution  is  find  a  process  that  prevents  excessive  mixing  without  increasing 
the  anisotropy  of  the  TKE  budget. 
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COMPARISON  OF  MODELS  AT  OWS  P.  19*, 
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COMPARISON  MODELS  AT  OWS  P,  1961 
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COMPARISON  OF  MODELS  AT  OWS  N.  1961 
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V.  SUMMARY  AND  RECOMMENDATIONS 


This  work  introduced  rotation  stress,  a  mechanism  by  which  planetary  rotation 
and  wind  direction  affect  geophysical  turbulent  boundary  layers.  Rotation  stress  is  the 
interaction  of  the  northward  component  of  planetary  rotation  and  the  east-west 
Reynolds  stress.  This  interaction  causes  a  redistribution  of  turbulent  kinetic  energy 
(TKE)  from  the  east-west  to  the  vertical  component  for  easterly  winds,  or  from  the 
vertical  to  the  east-west  component  for  westerly  winds.  This  exchange  of  TKE 
between  the  components  results  from  inertial  motions  of  the  turbulence  in  the  x-z 
plane. 

The  effects  of  rotation  stress  were  compared  with  those  that  result  from  the 
rotational  enhancement  to  the  TKE  dissipation  that  was  proposed  by  Garwood  (1977). 
There  are  two  major  results.  First,  the  TKE  distribution  is  more  isotropic  and  in  better 
agreement  with  laboratory  results  when  rotation  stress  is  included  in  a  TKE  model 
than  for  the  model  that  includes  rotationally  enhanced  dissipation  or  for  the  model 
with  no  planetary  rotation  effects.  Secondly,  rotation  stress  alters  the  frequency  and 
magnitude  of  diurnal  to  synoptic  scale  mixing  events  during  late  winter  through  early 
spring  at  mid  latitudes. 

The  ratio  of  the  vertical  to  the  total  TKE  is  a  measure  of  the  isotropy  of  the 
TKE  distribution  between  the  vertical  and  the  horizontal  components.  The  ratio 
would  be  0.33  if  the  TKE  distribution  was  isotropic.  Laboratory  measurements  of 
wall-bounded,  neutrally  stratified  shear  flows  indicate  that  the  ratio  is  0.16  to  0.3.  For 
stably  stratified  shear  flows  such  as  the  upper  ocean  mixed  layer,  the  ratio  would  be 
smaller.  For  late  winter  to  early  spring  at  OWS  P,  which  is  predominately  forced  by- 
surface  shear  stress,  typical  values  of  the  ratio  of  vertical  to  total  TKE  are  less  than 
0.03  for  the  rotationally  enhanced  dissipation  model  compared  to  an  average  of  0.1  for 
the  rotation  stress  model.  The  nondimensional  studies  show  that  this  ratio  can  be 
more  sanable  as  a  function  of  rotation  stress  than  as  a  function  of  rotationally 
enhanced  dissipation  since  rotation  stress  can  cause  a  transfer  of  energy  between  the 
her,  r'  .I  and  \ crtical  TKE  components. 

Him  ';ie  equilibrium  depth  and  the  entrainment  rate  are  functions  of  rotation 
■  f —  a  -Tt.c  stress,  diurnal  ^hallowing  events  occur  more  frequently  for 
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westerly  wind  regimes  and  the  equilibrium  depth  s  in  :o  40'  .  *nalle*er  tV  •  ••  •  •  .c 
model  which  neglects  planetar,  rotation  On  svnoptii.  Mme  vcjIcv  the  m.'.v 
simulated  with  the  rotation  stress  modei  typically  deepens  |o  .  mere  (<:r  ea^c  •  .< , -u. * 
or  about  10%  to  15“  >  less  for  westerly  winds  than  it  is  lor  the  mi -del  th.it  neglect 
planetary  rotation. 

For  late  winter  to  early  spring,  typical  values  ol  ’he  entrainment  rate  !roin  the 
simulations  with  rotation  stress  were  twice  as  large  at  OWS  \  and  lour  times  as  arge 
at  OWS  P  as  for  the  simulations  with  the  rotationjilv  enhanced  dissipation  V'  j 
result,  annual  simulations  using  the  rotationally  enhanced  JisMp.r-.cn  nad  sha!!i  s-er 
average  winter  mixed  layer  depths  than  those  for  rotation  stress  ihc»e  sh.ido'cr 
winter  mixed  layer  depths  better  fit  the  seasonal  trend  ol  the  OWS  ^h*cr\jti<  ns 
However,  the  functional  form  of  the  entrainment  rate  as  a  function  of  rotation  s’ress  .* 
consistent  with  the  need  for  decreased  mixing  efficiency  for  deep  mixed  .avers  and  large 
values  of  the  surface  stress.  This  is  not  the  wase  for  the  rotationa.  enhanced  dissipation 
rate.  The  entrainment  rate  decreases  to  zero  with  increasing  rotation  stress,  u.-.ereas  n 
approaches  an  asymptotic  value  greater  than  zero  as  the  rotationally  enhancement  to 
the  dissipation  increases.  Therefore,  the  rotation  stress  could  prevent  mixed  ’aver 
deepening  for  sufficiently  large  values  of  westerly  wind  or  deep  mixeJ  ;av  ers 
Conversely,  the  entrainment  rate  never  goes  to  zero  for  the  rotational:-  enhanced 
dissipation  case. 

At  midlatitudes,  the  rotation  stress  does  not  consistently  prevent  cxcessivclv  Jeep 
wintertime  mixed  layers.  Some  additional  mechanism  is  needed  to  reduce  mixing 
efficiency  for  deep  mixed  layers  and,  or  large  surface  stresses  at  mid  latitudes 
Enhanced  dissipation  does  prevent  excessively  deep  wintertime  mixed  layer  depths; 
however,  the  TKE  distribution  then  becomes  very  anisotropic.  The  degree  of 
anisotropy  appears  to  be  excessive  based  on  laboratory  measurements,  but  there  are  no 
oceanic  measurements  to  confirm  this. 

The  temporal  and  spatial  resolution  and  the  precision  of  the  OWS  data  sets  used 
in  this  study  were  not  sufficient  to  confirm  the  increased  frequency  and  magnitude  of 
diurnal  to  synoptic  scale  shallowing  and  deepening  events  predicted  by  the  rotation 
stress  model.  Data  from  experiments  such  as  the  Storm  Transfer  and  Response 
Experiment  (STREX)  and  the  Mixed  Layer  Dynamics  Experiment  (M1LDEX)  may 
be  sufficiently  accurate  with  enough  temporal  resolution  to  test  the  effects  of  rotation 
stress.  However,  these  data  sets  may  not  contain  the  contrasting  cast-west  wind  stress 
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events.  which  is  a  requirement  to  te«*  the  rotation  stress  diet's  It  was 'm  pcJ  that  the 
OWS  station  data  sett  would  contain  sufficient  events  to  -.onf;rm  the  effects  ul  rot  at  .on 
tt'est  Alto  the  O'*  S  data  provided  the  opportune  to  studs  the  seasons!  variability 
of  the  effects  of  planetary  rotation  on  the  upper  ocean  The  study  docs  indicate  that 
the  planetars  rotation  ha»  the  greatest  effect*  for  winter  and  earls  spring  m:\cd  lasers 
I  he  impact  of  planetary  rota\on  :s  the  smallest  lor  summer  mixed  lasers 

The  greatest  differences  among  the  simulations  :n  this  studs  appear  in  the  I  K I 
budgets  calculated  bs  the  models  !he  sanations  m  the  mixed  laser  temperature  are 
relum el\  small  due  to  the  large  heat  cjpacits  ol  water  Hence,  the  ddTerer.ccs  among 
he  predicted  and  observed  temperatures  are  small  The  differences  in  the  mixed  laser 
depths  predicted  by  the  models  are  more  obvious.  However,  the  models  -an  Nr  tuned 
to  predict  sin  .sod  laser  depths  even  though  the  I  K 1  budgets  of  the  models  are 
quite  different  »iso  the  diffci cnees  between  the  turbulent  boundary  layer  depth 
predicted  by  the  models  and  the  depth  of  the  well-mixed  laser  defined  by  the 
temperature  observations  obscure  the  comparison  between  the  simulations  and  the 
observations.  The  turbulent  boundary  layer  depth  and  the  depth  of  the  isothermal 
layer  in  the  upper  ocean  do  not  always  coincide,  and  the  OWS  data  are  not  alwass 
sufficient  to  define  the  isothermal  layer  accurately  or  consistently 

This  studs  has  shown  that  substantial  variations  exist  in  the  isotropy  of  the  !  Kl 
distribution  predicted  by  different  TKF.  models  of  the  upper  ocean  The  results  ol  this 
studs  suggest  that  parameterization*  of  various  processes  are  not  complete.  In 
particular,  the  parameterization*  for  dissipation,  shear  stress  and  the  prcsvurc-ratc  of 
strain  require  further  theoretical  and  observational  work.  The  importance  of  rotation 
stress  has  been  indicated  by  this  study,  but  the  data  were  insufficient  to  verily 
completely  the  effects  of  rotation  stress. 

The  results  of  this  research  focus  attention  on  specific  observational  requirements 
that  need  to  be  met  in  order  to  properly  evaluate  and  improve  turbulence  models  lor 
applications  to  the  oceanic  environment.  Direct  measurements  of  TKE  budget  in  the 
upper  ocean  are  required  to  test  the  current  hypotheses  of  ocean  planetary  boundary 
layer  dynamics.  A  field  experiment  that  focuses  on  the  upper  ocean  turbulence  from 
the  sea  surface  to  the  transient  thermocline  is  necessary  to  provide  this  data.  The  first 
phase  cf  such  an  idealized  program  should  concentrate  on  vertical  exchanges  and 
budgets  in  an  open  ocean  area  with  relatively  small  but  well-observed  horizontal 
variability.  Measurements  must  include  the  surface  heat,  salt  and  vector  momentum 
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(luxes.  Estimates  from  bulk  formulas  would  probably  not  be  sufficiently  accurate.  The 
mean  temperature,  salinity  and  velocity  profiles  from  the  surface  through  the  transient 
thcrmocline  must  also  be  measured  simultaneously.  The  vertical  and  temporal 
variations  and  the  vertical  shear  of  the  normal  Reynolds  stresses  (TKE  components) 
are  the  turbulence  quantities  to  be  measured  in  the  first  phase.  The  vertical  shear  of 
the  horizontal  and  vertical  turbulent  velocities  have  previously  been  measured  to 
estimate  dissipation.  The  second  phase  would  extend  the  first  phase  measurements  to 
include  the  vertical  and  temporal  vaiations  of  the  shearing  Reynolds  stresses  and  the 
vertical  buoyancy  flux  through  the  mixed  layer  and  the  entrainment  zone.  The  third 
phases  would  begin  to  incorporate  three-dimensional  effects  and  the  interactions 
between  the  turbulence  and  aspects  of  the  mean  (low  such  as  upwelling  and  internal 
waves. 

Practical  forecasts  of  mixed  layer  temperature  and  depth  can  not  be  improved 
until  the  turbulence  in  the  OPBL  is  better  understood  and  better  parameterized  in  the 
models  of  the  OPBL.  This  impacts  both  short-term  forecasts  of  the  upper  ocean  that 
are  important  for  improved  acoustical  forecasts  and  weather  forecasts  required  for 
naval  operations  and  climate  forecasts  of  the  ocean-atmosphere  system. 


APPENDIX 

FIGURES 


Vertical  North 


Elongation 


Figure  1.1.  The  vertical  elongation  or  contraction  of  turbulent  eddies  in  the 
castward-verticai  plane  by  the  northward  component  of  planetar*'  vorticitv. 
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Figure  4.1c.  The  distribution  of  XBT  casts  at  OWS  P  for  1966.  The  class  interval  is 
22.8  davs. 
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Figure  4.  Id.  The  distribution  of  XBT  casts  at  OWS  P  for  1967.  The  class  interval  is 
22. S  days. 
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Ficure  4.1e.  The  distribution  of  XBT  casts  at  OWS  N  for  1961.  The  class  interval  is 
22*3  days. 
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figure  4.3.  cont.  g)  total  1  Kl:  (curs'2),  h)  ratio  of  vertical  to  total  TKH  and  i)  TKE 
transfer  rate  due  to  the  pressure-rate  of  strain  minus  rotation  stress  (enrs'" ). 
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f  igure^4.6.  1967  Days  250-257,  OWS  P.  a)  mixed  layer  depth  (m),  b)  dissipation  rate 
(enrs'- )  and  c)  entrainment  rate  (enrs*3)  (or  the  RSTAR  (solid),  IIS  TAR  (dotted)  and 
ZSI AR  (dashed)  models.  Ihe  X  s  are  the  mixed  layer  depths  calculated  from  the  BT  s. 
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Figure  4.8a.  Annua!  temperature  bias  values  during  the  tuning  of  the  RSTAR  model 
with  p2«0.5  for  the  entire  year  of  1965  at  OWS  P. 
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fiL’iirc  4.10b.  OWS  P.  1961.  Mixed  layer  temperature:  top  curve  is  the  observations, 
middle  curve  is  the  RSTAR  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  4.11b.  OWS  P,  1965.  Mixed  layer  temperature:  top  curve  ts  the  observations, 
middle  curve  is  the  RSI  AB.  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  4.12b.  OWS  P.  1966,  Mixed  layer  temperature:  tep  curve  is  the  observations, 
middle  curve  is  the  RSTAR  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  4.13a.  OVVS  P.  1967,  Mixed  layer  depth:  top  curve  is  the  observations,  middle 
curve  is  the  RSTAR  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  ^.13b.  OWS  P.  1967.  Mixed  layer  temperature:  top  curve  is  the  observations. 
iniJcie  curve  is  the  RSTAR  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  4.14a.  OWS  N,  1961.  Mixed  layer  depth:  top  curve  is  the  observations,  middle 
curve  is  the  RSTAR  model,  bottom  curve  is  the  ZSTAR  model. 
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Figure  4.14b.  OWS  X.  1961,  Mixed  layer  temperature:  top  curve  is  the  observations 
:udd;e  curve  is  the  RSTAR  modei,  bottom  curve  is  the  ZSTAR  model. 
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